- Email: [email protected]
Probabilistic robust parity relation for fault detection using polynomial chaos
Proceedings of 20th The International Federation of Congress Automatic Control Proceedings of the the 20th World World Congress Proceedings of the 20th9-14, World The International Federation of Congress Automatic Control Toulouse, France, July 2017 The International Federation of Automatic Control Available online at www.sciencedirect.com The International of Automatic Control Toulouse, France, July Toulouse, France,Federation July 9-14, 9-14, 2017 2017 Toulouse, France, July 9-14, 2017
ScienceDirect
IFAC PapersOnLine 50-1 (2017) 1019–1024 Probabilistic robust parity relation for fault Probabilistic robust parity relation for fault Probabilistic relation for fault detectionrobust using parity polynomial chaos detection using polynomial chaos detection using polynomial chaos ∗ ∗ ∗
Yiming Wan ∗ Eranda Harinath ∗ Richard D. Braatz ∗ Yiming Yiming Wan Wan ∗∗ Eranda Eranda Harinath Harinath ∗∗ Richard Richard D. D. Braatz Braatz ∗∗ ∗ Yiming Wan Eranda Harinath Richard D. Braatz Massachusetts Institute of Technology, Cambridge, MA, USA ∗ ∗ Massachusetts Institute of Technology, Cambridge, MA, USA ∗ Massachusetts Institute of Technology, Cambridge, MA, USA Massachusetts Institute Technology, Cambridge, MA, USA e-mail: {ywan,oferanda, braatz}@mit.edu e-mail: e-mail: {ywan, {ywan, eranda, eranda, braatz}@mit.edu braatz}@mit.edu e-mail: {ywan, eranda, braatz}@mit.edu Abstract: In this paper, a robust parity relation approach is proposed for fault detection Abstract: In paper, robust relation is for detection Abstract: In this thiswith paper, aprobabilistic robust parity parity relation approach approach is proposed proposed for fault fault detection of linear systems the a time-invariant parameter uncertainties. To deal with Abstract: In thiswith paper, aprobabilistic robust parity relation approach is proposed for fault detection of linear systems the time-invariant parameter uncertainties. To deal of linear systems with the time-invariant deal with with polynomial dependence on probabilistic uncertain parameters, a setparameter of parityuncertainties. relations withTopolynomial of linear systems with the probabilistic time-invariant parameter uncertainties. To deal with polynomial dependence on uncertain aa set of parity polynomial dependence on to uncertain parameters, set of condition. parity relations relations with polynomial parameterization is derived decoupleparameters, the unknown initial Due towith the polynomial polynomial dependence on to uncertain parameters, a set of condition. parity relations with polynomial parameterization is derived decouple the unknown initial Due to the parameterization derived to decouple the unknown Due to the polynomial polynomial structure having aisfinite degree, the generated residual initial vector condition. admits an exact polynomial chaos parameterization isfinite derived to decouple the unknown initial condition. Due to the polynomial structure having a degree, the generated residual vector admits an exact polynomial chaos structure a finite generated vectoronline admits an exact polynomial chaos expansionhaving with the samedegree, degree,the which enablesresidual an efficient computation for quantifying structure having a finite degree, the generated residual vector admits an exact polynomial chaos expansion with degree, which an efficient computation for expansion with the same degree, which enables an moments efficient online online computation foraquantifying quantifying the moments of the the same residual vector. Theenables obtained are used to calculate confidence expansion with the same degree, which enables an efficient online computation for quantifying the moments of the residual The obtained moments used calculate aa confidence the moments of the the consistency residual vector. vector. Thethe obtained moments are used to to calculate confidence set for checking between observed system are behavior and the parity relations. the moments of the the residual vector. Thethe obtained moments are used to calculate a confidence set for between observed system and relations. set for checking checking the consistency consistency between therelation observed system behavior behavior and the the parity parity relations. Compared to deterministic robust parity methods, the employment of probabilistic set for checkingdeterministic the consistency between therelation observed system behavior and the parity relations. Compared robust the employment of Compared toallows deterministic robust parity parity methods, the alarm employment of probabilistic probabilistic informationto a less conservative way relation to ensuremethods, a low false rate while maintaining Compared to deterministic robust parity relation methods, the employment of probabilistic information a conservative to low rate maintaining information allows a less less conservative way to ensure ensure low false false alarm alarm rate while while maintaining a high fault allows detection rate. In contrastway to the existingaa polynomial chaos-based fault diagnosis information allows a less conservative way the to ensure a polynomial low false alarm rate while maintaining aaliterature, high detection In chaos-based fault high fault fault detection rate. In contrast contrast to the existing existing polynomial chaos-based fault diagnosis diagnosis our proposedrate. approach avoidsto truncation errors of polynomial chaos expansions, and aliterature, high fault detection rate. In contrast to the existing polynomial chaos-basedexpansions, fault diagnosis our approach errors polynomial and literature, our proposed proposed demanding. approach avoids avoids truncation errors of of is polynomial chaos expansions, and is less computationally The truncation proposed approach illustratedchaos by using a two-tank literature, our proposed demanding. approach avoids truncation errors of polynomial chaosusing expansions, and is less is less computationally computationally demanding. The The proposed proposed approach approach is is illustrated illustrated by by using aa two-tank two-tank system example. is less computationally demanding. The proposed approach is illustrated by using a two-tank system system example. example. system © 2017, example. IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Fault detection, parametric uncertainties, parity relation, polynomial chaos. Keywords: Keywords: Fault Fault detection, detection, parametric parametric uncertainties, uncertainties, parity parity relation, relation, polynomial polynomial chaos. chaos. Keywords: Fault detection, parametric uncertainties, parity relation, polynomial chaos. 1. INTRODUCTION have a low probability of occurrence. For tractable online 1. have of tractable 1. INTRODUCTION INTRODUCTION have aa low low probability probability of occurrence. occurrence. For tractable online computation, a convex set needs to be For adopted as anonline outer 1. INTRODUCTION have a low probability of occurrence. For tractable online computation, a convex set needs to be adopted as computation, a to convex set needs adopted as an an outer outer The parity relation approach, initiated by Chow and approximation the exact set to of be system behaviors precomputation, a to convex exact set needs to be adopted as an outer The parity relation approach, by set of behaviors preThe parity relation approach,asinitiated initiated by Chow Chow and approximation to the the exact setrelations, of system system behaviors preWillsky (1984), is recognized one of main classesand of approximation dicted by the uncertain parity which introduces The parity relation approach,asinitiated by Chow and approximation to the exact setrelations, of system behaviors preWillsky (1984), is one of uncertain parity which Willsky (1984), is recognized recognized asdetection one of of main main classes of dicted dicted by by the the that uncertain parity relations, which introduces introduces techniques for model-based fault (FD)classes (Gertler, conservatism further reduces FD performance. Willsky (1984),model-based is recognized asdetection one of main classes of dicted by the uncertain parity relations, which introduces techniques (FD) techniques for model-based fault detection (FD) (Gertler, (Gertler, conservatism that that further further reduces reduces FD FD performance. performance. 1998; Ding,for 2013). The parityfault relation is a receding horizon conservatism techniques for model-based fault detection (FD) (Gertler, conservatism that further of reduces FD performance. By explicit consideration the probabilistic 1998; Ding, 2013). The parity relation is a receding horizon 1998; Ding, 2013). parity relation isfrom a receding horizon By explicit consideration of the probabilistic information input-output (I/O)The model decoupled the unknown information 1998; Ding, 2013). The parity relation is a receding horizon By uncertainties, explicit consideration of the probabilistic information the probabilistic robust FD approach input-output (I/O) decoupled the unknown input-output (I/O) model decoupled from the or unknown initial states by usingmodel a so-called parityfrom matrix vector. of By explicit consideration of the probabilistic information of uncertainties, the probabilistic robust FD approach input-output (I/O) model decoupled from the or unknown of uncertainties, the probabilistic robust FD approach has received some attention for reducing the conservatism initial states by using a so-called parity matrix vector. initial by using so-called the parity matrix orbetween vector. of Faults states are detected bya checking consistency uncertainties, the probabilistic robust FD approach has received some attention for reducing the conservatism initial states by using a so-called parity matrix or vector. has received some attention for reducing the conservatism of the deterministic robust approaches (Esfahani and Faults are detected by checking the consistency between Faults are relation detectedand bythe checking the system consistency between has the parity measured I/O behavior. received some attention for reducing the conservatism the deterministic robust approaches (Esfahani and Faults are relation detectedand bythe checking the system consistency between of of the deterministic robust approaches (Esfahani and Lygeros, 2016; Zhong et al., 2016). As a computational tool the parity measured I/O behavior. the parity relationfor and theinmeasured system I/O behavior. A key challenge FD real applications is to ensure Lygeros, of the deterministic robust approaches (Esfahani tool and 2016; Zhong et al., 2016). As a computational the parity relation and the measured system I/O behavior. Lygeros, 2016; Zhong et al., 2016). As a computational tool for probabilistic uncertainty propagation, the generalized A key for FD to A key challenge challenge for uncertainties FD in in real real applications applications isalarm to ensure ensure robustness to model for low falseis rate Lygeros, 2016; Zhong et al., 2016). As a computational tool for probabilistic uncertainty propagation, the generalized A key challenge for uncertainties FD in real applications isalarm to ensure for probabilistic the employed generalized chaosuncertainty (PC) theorypropagation, (Xiu, 2010) was in robustness tomaintaining model for low false falserate rate polynomial robustness model uncertainties low alarm rate (FAR) whileto high faultfor detection (FDR). probabilistic uncertainty propagation, the employed generalized polynomial chaos (PC) theory (Xiu, 2010) was in robustness tomaintaining model uncertainties for low falserate alarm rate for polynomial chaos (PC) theory (Xiu, 2010) was employed in the probabilistic robust FD methods proposed by Mesbah (FAR) while high fault detection (FDR). (FAR) whiledevelopment maintaining of high detection rateapproach (FDR). polynomial The early thefault parity relation chaosrobust (PC) theory (Xiu, 2010) was employed in the probabilistic FD methods proposed by Mesbah (FAR) while maintaining high fault detection rate (FDR). the probabilistic robust FD methods proposed by Mesbah et al. (2014); Du et al. (2015, 2016). By expressing the The early development of the parity relation approach The early development of the parity relation approach mainly focused on additive unknown disturbances, and the probabilistic robust FD methods proposed by Mesbah et al. (2014); Du et al. (2015, 2016). By expressing the The early development of the parity relation approach et al. (2014); Du et al. function (2015, 2016). By expressing the dependence as a series of orthogonal mainly focused on unknown and mainly focused on additive additive unknown disturbances, disturbances, and uncertainty could not effectively address multiplicative uncertainties. al. (2014); Du et al. function (2015, 2016). By expressing the uncertainty dependence as aa series of orthogonal mainly focused on additive unknown disturbances, and et uncertainty dependence function as series of orthogonal polynomials, the PC framework captures the evolution of could not effectively address multiplicative uncertainties. could not effectively address multiplicative uncertainties. uncertainty dependence function as a series of orthogonal polynomials, the PC PC framework framework captures the evolution evolution of could not effectively address multiplicative uncertainties. the captures the of By employing set-based tools, the set-membership parity polynomials, probability distributions with the PC expansion (PCE) polynomials, the PC framework captures the evolution of By employing set-based tools, parity probability These distributions with the PC expansion expansion (PCE) By employing set-based tools, the the set-membership parity coefficients. probability distributions with PC (PCE) relation approach was proposed toset-membership deal with determinisPC-based FDthe methods rely on applying By employing set-based tools, the set-membership parity probability distributions with the PC expansion (PCE) relation approach was proposed to deal with determiniscoefficients. These PC-based FD methods rely on applying relation approach was proposed to deal with determiniscoefficients. PC-based FD methods rely on applying tic unknown-but-bounded parametric uncertainties (Ploix the intrusiveThese Galerkin projection to the differential equarelation approach was proposed to deal with determiniscoefficients. These PC-based FD methods rely on applying tic uncertainties (Ploix Galerkin projection to the equatic unknown-but-bounded parametric uncertainties (Ploix the the intrusive intrusive Galerkin projection to which the differential differential equaandunknown-but-bounded Adrot, 2006; Blesa et parametric al., 2012, 2016). The approach tions describing the plant dynamics, limits their aptic unknown-but-bounded parametric uncertainties (Ploix the intrusive Galerkin projection to which the differential equaand Adrot, 2006; Blesa et al., 2012, 2016). The approach tions describing the plant dynamics, limits their apand Adrot, 2006; Blesa et al., 2012, 2016). The approach tions describing the plant dynamics, which limits their apdecouples the unknown initial states by introducing the plicability to stable systems only, because the convergence and Adrot,the 2006; Blesa et al., 2012, 2016). The approach tions describing the systems plant dynamics, whichthe limits their apdecouples unknown initial by the stable because decouples the as unknown initial states by introducing introducing the plicability plicability todoes stable only, because the convergence convergence parity matrix the function of states the uncertain parameters. of the PCEto notsystems hold foronly, unstable systems. Moreover, decouples the as unknown initial states by introducing the plicability todoes stable systems only, because the convergence parity matrix function of uncertain parameters. PCE not for Moreover, parity matrix in as the the function of the the uncertain parameters. of of the the PCE does to notbehold hold for unstable unstable systems. Moreover, This results a set of parity relations parameterized the PCE needs truncated, and systems. no systematic way parity matrix in as the function of the uncertain parameters. of the PCE does to notbehold for unstable systems. Moreover, This results a set of parity relations parameterized the PCE needs truncated, and no systematic This in aparameters, set of parity relations theavailable PCE needs to be truncated, and no systematic way by theresults uncertain rather than aparameterized single parity is to account for such truncation errors in way the This results in aparameters, set of parity relations parameterized theavailable PCE needs to be truncated, and no systematic way by the rather than parity to for by the uncertain uncertain rather than aa single single parity is available to account account for such such truncation truncation errors errors in in the the relation as in the parameters, case of additive disturbances. Then the is robust FD system design. by the uncertain parameters, rather than a single parity availablesystem to account for such truncation errors in the relation as case of Then the is robust relation as in in isthe the case by of additive additive disturbances. the robust FD FD system design. design. FD decision made checkingdisturbances. whether the Then observed relation as in isthe case by of additive disturbances. Then the robustwork FD system design.the robust FD problem for lininvestigates FD decision whether the observed FD decision is made made by checking checking whether observed system behavior is admitted by the set ofthe parity rela- This Thissystems work investigates investigates the robust robust FD problem problem for linlinFD decision is made by checking whether the observed This work the FD for ear with polynomial dependence on probabilissystem behavior is admitted by the set of parity relasystem behavior such is admitted by the robust set of parity rela- This tions. Although a deterministic work investigates the robust FD problem for linear systems with polynomial dependence on probabilissystem behavior is admitted by the set of parity relaear systems with polynomial dependence on probabilistic time-invariant parametric uncertainties. The above tions. Although such aa deterministic robust parity tions. Althoughcan such deterministic robust even parityinrelarelation approach ensure zero false alarms the ear systems with polynomial dependence on The probabilistic parametric above tions. Although such a deterministic robust even parity relatic time-invariant time-invariant parametric uncertainties. uncertainties. The robust above of the deterministic and probabilistic tion ensure false tion approach can ensure zero false alarms alarms even in in the the limitations worstapproach case, the can resulting FDzero performance is compromised, tic time-invariant parametric uncertainties. The above limitations of the deterministic and probabilistic robust tion approach can ensure zero false alarms even in the FD limitations of the deterministic and probabilistic robust methods can be effectively handled by our proposed worst case, the resulting FD performance is compromised, worst case, the worst-case resulting FDparameter performance is compromised, because some uncertainties often limitations of the deterministic and probabilistic robust can effectively handled by proposed worst case, the worst-case resulting FDparameter performance is compromised, FD methods methods robust can be be parity effectively handled by our ourusing proposed probabilistic relation approach PC. because some uncertainties often FD because some worst-case parameter uncertainties often FD methods robust can be parity effectively handled by ourusing proposed ⋆ probabilistic relation approach PC. Support some from Novartis is acknowledged. because worst-case parameter uncertainties often probabilistic robust parity relation approach using PC. ⋆ probabilistic robust parity relation approach using PC. ⋆ Support Support from from Novartis Novartis is is acknowledged. acknowledged.
⋆ Support from Novartis is acknowledged. Copyright © 2017 IFAC 1042 2405-8963 © 2017, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Copyright © 2017 IFAC 1042 Copyright ©under 2017 responsibility IFAC 1042Control. Peer review of International Federation of Automatic Copyright © 2017 IFAC 1042 10.1016/j.ifacol.2017.08.210
Proceedings of the 20th IFAC World Congress 1020 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
Similarly to Ploix and Adrot (2006); Blesa et al. (2016), a parity matrix is designed as a polynomial function of the uncertain parameters. By exploiting the polynomial parametric dependence, both the uncertain parity relation and the generated residual vector are written into a linear form with respect to the monomials of uncertain parameters, which enables efficient online stochastic uncertainty quantification using offline PCE of the monomials. Moreover, due to the polynomial parametric dependence, the involved PCEs are exact without any truncation errors. Using the moment information computed from the PCE of the residual vector, a Chebyshev bound is calculated for checking the consistency between the parity relation and the observed I/O data. This paper is organized as follows. Section 2 states the probabilistic FD problem. Section 3 presents the computation of the polynomial parity matrix. The residual uncertainty quantification is discussed in Section 4. Section 5 constructs the consistency test for the FD decision making. The simulation results are presented in Section 6. 2. PROBLEM STATEMENT 2.1 System description Consider the linear uncertain dynamical system x(k + 1) = A(θ)x(k) + B(θ)u(k) y(k) = C(θ)x(k) + D(θ)u(k) + e(k)
(1)
where the system inputs u(k) ∈ Rnu and outputs y(k) ∈ Rny are measurable, the states x(k) ∈ Rnx are unmeasured, e(k) denotes the additive stochastic measurement noises with zero mean and known covariance Σe , and A(θ), B(θ), C(θ), and D(θ) are matrices whose entries are polynomials of the uncertain parameters θ ∈ Rnθ . Such parametric uncertainties are often encountered in practice when specific knowledge is available on how uncertain parameters enter a system model (Ackermann, 2002). The uncertain parameters θ are time invariant, which is also a good approximation for slowly time-varying parameters (Blesa et al., 2016). We also assume that the entries of θ are mutually independent random variables with known time-invariant probability density functions (PDFs) and supports. Such a description of θ can be derived from the a prior knowledge of the underlying system, or offline parameter identification from data.
yk,h
y(k − h + 1) y(k − h + 2) , = .. .
C(θ) C(θ)A(θ) Oh (θ) = .. .
,
C(θ)Ah−1 (θ) 0 ··· 0 . .. . .. C(θ)B(θ) D(θ) Hh (θ) = .. .. .. . 0 . . C(θ)Ah−2 (θ)B(θ) C(θ)Ah−3 (θ)B(θ) · · · D(θ) (3) and uk,h and ek,h are defined similarly to yk,h . A parity matrix V (θ) is designed to satisfy V (θ)Oh (θ) = 0 for any θ (4) so that the premultiplication of (2) by this matrix generates the parity relation (5) V (θ)(yk,h − Hh (θ)uk,h − ek,h ) = 0, without the unknown initial term Oh (θ)x(k − h + 1). The parity matrix V (θ) is assumed to be a polynomial matrix function of the uncertain parameters θ, due to the polynomial dependence of Oh (θ) on θ. Using this parity relation (5), the residual vector is generated with polynomial dependence on θ and e(k): r(k, θ, e) = G(θ)zk,h − V (θ)ek,h (6) [ ⊤ ] ⊤ ⊤ with G(θ) = [V (θ) −V (θ)Hh (θ)] and zk,h = yk,h uk,h .
y(k) D(θ)
With the bounded descriptions θ ∈ Θ and e(k) ∈ E, the set-based consistency test was proposed as (Ploix and Adrot, 2006; Blesa et al., 2016) (7) 0 ∈ Γk , where Γk = {r|r = G(θ)zk,h − V (θ)ek,h , θ ∈ Θ, e(k) ∈ E}. When (7) holds, there exists θ0 ∈ Θ and {e0 (i) ∈ E}ki=k−h+1 such that r(k, θ0 , e0 ) = 0 and the occurrence of faults cannot be concluded. Otherwise, the observed I/O data are inconsistent with the parity relation, which indicates the presence of faults. The exact computation of the set Γk in (7) is often prohibitive due to its nonlinear dependence on θ. This difficulty was circumvented in Blesa et al. (2016) by adopting a different linear parametrization, at the cost of reduced FDR. Since the residual r(k, θ, e) is a polynomial vector of θ and e, (6) was equivalently rewritten as (Blesa et al., 2016) r (k, ξ(θ), n(θ, ek,h )) = η(k, ξ(θ)) + n(θ, ek,h ), (8) ⊤
where ξ(θ) = [ξ0 (θ) · · · ξNG −1 (θ)] , η(k, ξ(θ)) =
2.2 Robust parity relation
N∑ G −1
φi (zk,h )ξi (θ),
(9)
i=0
With the unknown-but-bounded assumption on the parametric uncertainties θ and the measurement noises e(k) in the system (1), Ploix and Adrot (2006) proposed a deterministic robust parity relation approach which is briefly reviewed below to provide background and to motivate our investigated problem. From the system model (1), the stacked output equation over a time window [k − h + 1, k] can be derived as where
yk,h = Oh (θ)x(k − h + 1) + Hh (θ)uk,h + ek,h
(2)
(10) n(θ, ek,h ) = −V (θ)ek,h , NG represents the number of monomials in the polynomial matrix G(θ), ξi (θ) denotes the monomial basis, and the corresponding monomial coefficient φi (zk,h ) is linear in the I/O data zk,h . With the new parameters ξ and n, the set Γk in (7) is replaced by (11) Υk = {r | r = η(k, ξ) + n, ξ ∈ Ξ, n ∈ N } . For the sake of tractable computation, the sets Ξ and N are used as outer approximations of {ξ|ξ = ξ(θ), θ ∈ Θ} and {n|n = V (θ)ek,h , θ ∈ Θ, e(k) ∈ E}. Therefore, Γk ⊆ Υk , and using 0 ∈ Υk as the consistency test reduces FDR.
1043
Proceedings of the 20th IFAC World Congress Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
Being a worst-case formulation, the above deterministic parity relation approach does not take into account the probability of occurrence for all possible values of parametric uncertainties and noises, even for values that have vanishingly low probability (e.g., could occur once in 10100 years). As a result, the deterministic robust approach is conservative in ensuring FAR, thus compromising the FDR. This article employs PC theory in the robust parity relation so that probability information of parametric uncertainties in the linear system (1) can be incorporated to enhance FDR while ensuring a low FAR. 3. COMPUTING A POLYNOMIAL PARITY MATRIX This section presents how to compute the polynomial parity matrix V (θ). For the ease of understanding, we discuss first the case of a univariate uncertain parameter θ ∈ R , and then proceed to the general case of the multivariate uncertain parameters θ ∈ Rnθ . For a scalar θ, assume that the degrees of the polynomial matrices A(θ) and C(θ) are dA and dC , respectively. Then the observability matrix Oh (θ) has a degree of dO = dC + (h − 1)dA . (12) When a polynomial parity matrix V (θ) is selected with a degree dV , (4) becomes ) d (d V O ∑ ∑ Vi θ i Oj θj = 0. (13) i=0
j=0
The multiplication on the left-hand side of (13) is a polynomial with a degree dV + dO + 1, whose coefficient ∑i of the monomial basis with a degree i is j=0 Vj Oi−j . All of these monomial coefficients need to be zero, i.e., ∑i j=0 Vj Oi−j = 0, so that (13) holds for all θ. Equation 13 can be compactly written as O 0 O1 O2 · · · 0 0 O0 O1 · · · 0 0 0 O ··· 0 = 0. 0 (14) [V0 V1 · · · VdV ] �� � � ... ... ... . . . ... V 0 0 0 · · · OdO � �� � M
The multiplication of V and the ith block-column of M is the coefficient of the monomial basis with a degree i − 1. The maximum number of block-columns in M is determined by the degree of V (θ)Oh (θ), and the size of M is (dV + 1)hny × nM , with nM ≤ (dV + dO + 1)nx . This idea can be extended to the case of a nθ -dimensional ⊤ parameter vector θ = [θ1 θ2 · · · θnθ ] . Let θ[S] represent ∑nθ snθ s1 s2 a monomial basis θ1 θ2 · · · θnθ with a degree i=1 si and S = (s1 , s2 , · · · , snθ ). Then a ∑ polynomial matrix M (θ) can be decomposed as M (θ) = M[Si ] θ[Si ] , where M[Si ] is Si ∈S
the coefficient matrix for the monomial basis θ[Si ] , and the set S is determined by all the monomials presented in the polynomial matrix M (θ). If the degree of M (θ) is dM , the maximum number of monomial terms in M (θ) is (dM + nθ )! . (15) g(dM , nθ ) = dM !nθ !
1021
Based on the definitions, (4) becomes ( ) ∑ ∑ V[Si ] θ[Si ] O[Qi ] θ[Qi ] = 0 Si ∈S
(16)
Qi ∈Q
for multivariate uncertain parameters. Finding a polynomial parity matrix V (θ) again relies on solving VM = 0, with V and M defined similarly to (14). Now the size of M is g(dV , nθ )hny × nM , with nM ≤ g(dV + dO , nθ )nx . The number of block-columns nM can be smaller than g(dV +dO , nθ )nx due to the possible zero coefficients O[Qi ] . The existence of a solution to VM = 0 is ensured by having the matrix M have more rows than columns, i.e., g(dV , nθ )hny > g(dV + dO , nθ )nx . Using (15), this inequality can be rewritten as ) nθ ( ∑ dO nx > 1− , (17) d + d + i hn V O y i=1
which implies that by selecting a sufficiently high degree dV and a large length h of the time window, a polynomial parity matrix V (θ) exists so that (14) holds. 4. QUANTIFYING RESIDUAL UNCERTAINTY With the parity matrix obtained in Section 3, the residual vector can be generated from (6). As a function of the random parameters θ and the stochastic noises e(k), the generated residual vector is also stochastic. In this section, PC theory is applied for efficient online uncertainty quantification of the residual vector. 4.1 PCE of multivariate polynomials Consider a multivariate polynomial matrix denoted by m ∑ M (θ) = Mi ξi (θ), (18) i=0
with ξi (θ) representing a monomial basis with its degree at most p. Then M (θ) has an exact PCE representation M (θ) =
Np ∑
Ψj Φj (θ),
(19)
j=0
where Φj (θ) denotes the multivariate PC basis function. The total number of terms in (19) is Np + 1 = (nnθθ+p)! !p! , depending on the dimension nθ of θ and the order p of M (θ) (Xiu, 2010). The PCE coefficients {Ψj } for M (θ) can be efficiently computed by following two steps: (i) determine the PCE coefficients {βij } of each monomial basis ξi (θ), i.e., ξi (θ) =
Np ∑
βij Φj (θ);
(20)
j=0
and (ii) compute the PCE coefficients {Ψj } in (19) as m ∑ Mi βij (21) Ψj = i=0
by substituting (20) into (18). Example 1. The Legendre polynomials form an orthogonal basis with respect to a uniform random scalar θ over
1044
Proceedings of the 20th IFAC World Congress 1022 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
[−1, 1]. The first four Legendre polynomials are given by (Maˆıtre and Knio, 2010) 1 1 0 0 0 Φ0 (θ) Φ1 (θ) 0 1 0 0 θ Φ (θ) = − 1 0 3 0 θ2 . 2 2 2 Φ3 (θ) 0 − 32 0 52 θ3
By inverting the coefficient matrix, the monomials θ, θ2 , and θ3 can be represented as Φ0 (θ) θ 0 1 0 0 Φ1 (θ) θ 2 = 1 0 2 0 Φ (θ) , 3 3 2 3 3 2 0 5 0 5 θ Φ (θ) 3
where the elements of each row of the coefficient matrix are the PCE coefficients of the corresponding monomial. Then the PCE coefficients of any polynomial matrix (18) with a degree less than 4 can be determined according to (21).
The PC basis functions belong to the Askey scheme of orthogonal polynomials, satisfying (Xiu, 2010) ⟨Φi (θ), Φj (θ)⟩ = ⟨Φ2i (θ)⟩δij , where δij is the Kronecker delta, and ⟨·, ·⟩ defines the inner product ∫ ⟨α(θ), β(θ)⟩ = α(θ)β(θ)fθ dθ
respectively. The polynomial dependence of η(k, ξ(θ)) and n(θ, ek,h ) on θ in (24) allows the uncertainty quantification to be decomposed into four sequential steps: (i) For the selected polynomial bases {Φj (θ)}, compute the inner products {⟨Φ2j (θ)⟩} and the PCE coefficients of each monomial basis ξi (θ) as in (20). (ii) Compute Cov{n(θ, ek,h )}. From (10), it follows that ⊤ Cov{n(θ, ek,h )} = E{V (θ)ek,h e⊤ k,h V (θ)}
(25) = E{V (θ)Σe,h V ⊤ (θ)} = E{Σn (θ)}, where Σe,h = Ih ⊗ Σe denotes the covariance of ek,h . Note that Σn (θ) = V (θ)Σe,h V ⊤ (θ) in (25) has polynomial dependence on θ. After calculating the PCE coefficients of Σn (θ) from (21), E{Σn (θ)} is determined by its zero-order PCE coefficient according to (22). (iii) Compute E{η(k, ξ(θ))} and Cov{η(k, ξ(θ))} by using the online data zk,h . For a given time window, the I/O data zk,h can be regarded as deterministic. Then, by substituting (20) into (9), the PCE of η(k, ξ(θ)) is η(k, ξ(θ)) =
E{M (θ)} =
Np ∑ j=0
Cov{M (θ)} =
Ψj ⟨Φj (θ), 1⟩ = Ψ0 ,
N p Np ∑ ∑ i=1 j=1
=
Np ∑ j=1
(22)
Ψi Ψ⊤ j ⟨Φi (θ), Φj (θ)⟩
2 Ψj Ψ⊤ j ⟨Φj (θ)⟩.
(23)
Note that the inner products of basis functions {⟨Φ2j (θ)⟩} 1 for an uniform can be pre-computed, e.g., ⟨Φ2j (θ)⟩ = 2j+1 random scalar θ, so that the covariance in (23) can be efficiently calculated using PCE coefficients. 4.2 Polynomial chaos approach to computing residual mean and covariance In this subsection, the mean and covariance of the residual vector (6) is computed by applying PC theory reviewed in Section 4.1. For this purpose, the linear parameterization form (8), instead of (6), is adopted. Due to the probabilistic independence between θ and e(k), E{n(θ, ek,h )} = 0 and E{n(θ, ek,h )η ⊤ (k, ξ(θ))} = 0 are obtained, thus the mean and covariance of the uncertain residual vector (8) are µk = E{η(k, ξ(θ))}, and (24) Σk = Cov{η(k, ξ(θ))} + Cov{n(θ, ek,h )},
ηj (k)Φj (θ),
j=0
Ω
with respect to the PDF fθ of θ over the support Ω. The selection of PC bases is determined by probability distributions of random parameters, e.g., Legendre polynomials for uniform distributions, Hermite polynomials for Gaussian distributions. Exploiting the above orthogonality, the mean and covariance of the polynomial matrix M (θ) can be conveniently determined from its PCE coefficients {Ψj } by (Xiu, 2010)
Nq ∑
ηj (k) =
N∑ G −1
(26) φi (zk,h )βij ,
i=0
with q representing the degree of η(k, ξ(θ)). The PCE coefficients ηj (k) can be efficiently computed by using the online data zk,h and the PC coefficients βij obtained offline in (20). From the PCE coefficients ηj (k) in (26), E{η(k, ξ(θ))} and Cov{η(k, ξ(θ))} can be calculated from (22) and (23). (iv) According to (24), compute the mean µk and covariance Σk of the residual vector (8) using results from Steps (ii) and (iii). The first two steps can be performed offline, and the last two steps are online computations involving only basic matrix manipulations, which results in much more efficient computation, compared to using Galerkin projection based methods in Mesbah et al. (2014); Du et al. (2015, 2016). The involved PCEs are all exact with no truncation errors, due to the polynomial parametric dependence of a finite order in the system model (1). In contrast, the PCE-based methods in Mesbah et al. (2014); Du et al. (2015, 2016) have truncation errors that are not explicitly considered. 5. CONSISTENCY TEST In this section, the Chebyshev bound of the residual vector are calculated for the consistency test by using the residual mean µk and covariance Σk obtained in Section 4.2. Lemma 1. (Chen and Zhou (1997); Budny (2016)). For the stochastic residual vector rk = r(k, θ, e) with the mean µk and the covariance Σk ,
1045
(i) Assume Σk is positive definite with rank(Σk ) = nr . For any ϵ > 0, it holds that { } nr ⊤ , Pr (rk − µk ) Σ−1 k (rk − µk ) < ϵ > 1 − ϵ
Proceedings of the 20th IFAC World Congress Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
(ii) Assume Σk is rank deficient with rank(Σk ) = mk , and has a eigenvalue decomposition Σk = Uk Πk Uk⊤ , where the diagonal matrix Πk = diag(λ1 , · · · , λmk , 0, · · · , 0), with λi > 0. For any ϵ > 0, it holds that } { mk ⊤ , Pr (rk − µk ) Uk Tk Uk⊤ (rk − µk ) < ϵ > 1 − ϵ −1 where Tk = diag(λ−1 1 , · · · , λmk , 0, · · · , 0). For a nonsingular or singular covariance matrix Σk , Lemma 1 provides a confidence ellipsoid ∆k with a confidence level 1 − γ, i.e., Pr {rk ∈ ∆k } > 1 − γ, (27) where ∆k and γ are defined by { � } nr � ⊤ , or ∆k = rk �(rk − µk ) Σ−1 k (rk − µk ) < ϵ , γ = { � } ϵ m � k ⊤ , ∆k = rk �(rk − µk ) Uk Tk Uk⊤ (rk − µk ) < ϵ , γ = ϵ respectively. With a large ϵ, rk ∈ ∆k holds with a high confidence level 1 − γ, then the ellipsoid ∆k describes the most plausible region in which the residual vector lies. For this reason, the confidence ellipsoid ∆k can be used in the consistency test as Γk in (7), i.e., performing the consistency test by checking whether 0 ∈ ∆k holds. This approach is as simple as checking whether J < 1, with ⊤
J = ϵ−1 (rk − µk ) Σ−1 k (rk − µk ) , or
(28) ⊤ J = ϵ−1 (rk − µk ) Uk Tk Uk⊤ (rk − µk ) . The invalidation of J < 1 indicates the presence of faults. Remark 1. Since Lemma 1 holds for any PDF with the given mean and covariance, the confidence level 1 − γ in (27) is always conservative, i.e., the actual probability of rk ∈ ∆k can be much larger than 1 − γ. Therefore, it is possible to increase γ for an improved FDR without reducing FAR. 6. SIMULATION EXAMPLE Consider a two-tank process described by (Blesa et al., 2016) dx1 α1 √ α2 √ κ1 =− 2gx1 + 2gx2 + v1 dt A1 A1 A1 (29) dx2 α2 √ κ2 =− 2gx2 + v2 dt A2 A2 where α1 = α2 = 7 cm2 are the cross-sections of the two outlet holes, A1 = A2 = 28 cm2 are the cross-sections of the two tanks, κ1 = 2.33 cm3 /V s and κ2 = 1.34 cm3 /V s are the coefficients of the two inlet flows, v1 and v2 are input voltages to the pumps, x1 and x2 are the water levels of the two tanks, and g = 981 cm/s2 . The sampled ⊤ system inputs and outputs are u(k) = [v1 (k) v2 (k)] and ⊤ y(k) = [x1 (k) x2 (k)] + e(k), with the sampling interval ts = 1 s and the zero-mean white noises whose covariance is Σe = diag(10−4 , 10−4 ) cm2 . The water levels may slowly vary within the interval [xci − ∆xi , xci + ∆xi ], with xc1 = 30 cm, xc2 = 15 cm, and ∆x1 = ∆x2 = 5 cm. The continuous-time system (29) is discretized by the firstorder Euler method and linearized around the equilibrium point described by x01 , x02 , v10 , v20 , which results in the discrete-time model,
[
1023
] [ ] ][ δx1 (k + 1) a (x0 (k)) a12 (x02 (k)) δx1 (k) = 11 1 δx2 (k + 1) 0 a22 (x02 (k)) δx2 (k) ][ ] [ (30) δv1 (k) b11 0 , + 0 b22 δv2 (k)
where δvi = vi −vi0 , δxi = xi −x0i , and aij and bij are determined by constants of the system (29) and the uncertain water levels. To apply the proposed method, the nonlinear functions a11 (x01 (k)), a12 (x02 (k)), and a22 (x02 (k)) over the water level variation intervals [xci − ∆xi , xci + ∆xi ] were approximated by their Taylor expansions at the centers xc1 and xc2 , i.e., a ˆ11 (∆x1 (k)), a ˆ12 (∆x2 (k)), and a ˆ22 (∆x2 (k)), respectively. For this system, second-order Taylor expansions are sufficiently accurate. By doing so, the model (30) has polynomial dependence on two uncertain parameters θ1 = ∆x1 and θ2 = ∆x2 . Due to lack of in-depth probabilistic information, θ1 and θ2 are both assumed to be uniformly distributed in [−5, 5] cm, thus Legendre polynomial bases are adopted. The time window length h of the parity relation is selected to be h = 2, and a 2-by-4 polynomial parity matrix with a degree 2 is found. For both fault-free and faulty scenarios, 100 Monte Carlo (MC) runs are implemented with different time-invariant parameters θ in each run, with the fault being a leak in tank 1 injected at 20 s by adding ∆α1 = 0.01 cm2 to α1 in (29). Confidence ellipsoids ∆k and residual vectors of 100 MC runs with this fault are visualized in Fig. 1 for a sequence of 4 sampling instants, with γ = 0.15. The consistency condition 0 ∈ ∆k becomes invalid after 21 s, which is prompt detection of the leakage fault. The ensemble of the residual evaluations (28) in the fault-free and fault MC runs for the same γ shows zero FAR and high FDR accounting for the variations of uncertain parameters and measurement noises (Fig. 2). With a smaller leakage fault ∆α1 = 0.008 cm2 , the case of θ1 and θ2 being Gaussian distributed is also simulated. Each uncertain parameter has its standard deviations equal to 1 cm, so that its value almost surely belongs to the 5-sigma interval [−5, 5] cm. Two FD designs are considered: (i) uniform distribution over [−5, 5] cm is used for θ1 and θ2 in the design; and (ii) zero-mean Gaussian distribution with its standard deviation being equal to 1 cm is used in the design. For statistical performance evaluation, the FAR and FDR over a finite time window [k1 , k2 ] was used as defined by ∑k 2
number of false alarms at time instant i (number of fault-free MC runs) × (k2 − k1 + 1) ∑k2 i=k1 number of FDs at time instant i . FDR = (number of fault MC runs) × (k2 − k1 + 1) FAR =
i=k1
As shown in Figure 3, increasing γ from 0.1 to 0.26 leads to an improved FDR without sacrificing FAR, due to the reason explained in Remark 1. Moreover, the design (ii) has much higher FDR than the design (i) for all values of γ ∈ [0.1, 0.36], while it only suffers from higher FAR for γ = 0.3. This improvement of FD performance in the design (ii) benefits from exploiting the more informative probabilistic distribution information on the uncertain parameters.
1046
Proceedings of the 20th IFAC World Congress 1024 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
t = 20s residuals
0
-0.05 -0.1
confidence ellipsoids
0
r(2) k
r(2) k
0.05
REFERENCES
-0.05 -0.05
0
0.05
-0.05
0
r(1) k
t = 23s
0
r(2) k
r(2) k
0.05
r(1) k
t = 22s
0
Future work will focus on comparisons with existing setmembership parity relation methods, and extensions to nonlinear dynamical systems.
t = 21s 0.05
-0.05
-0.05 -0.1
-0.1 -0.05
0
0.05
-0.05
0
r(1) k
0.05
r(1) k
1
Normalized evaluation 0.5
Consistency evaluation
Consistency evaluation
Fig. 1. Confidence ellipsoids and residual vectors of 100 fault Monte Carlo runs from 20 to 23 s, with γ = 0.15.
0 10
15
6
Normalized threshold
20
25
30
35
40
45
50
45
50
Fault Injection
4 2 0 10
15
20
25
30
Time (s)
35
40
Fig. 2. The ensemble of the residual evaluations (28) in the fault-free and fault Monte Carlo runs for γ = 0.15.
FAR(%)
15 10
Uniform distribution Gaussian distribution
5 0 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34 0.36
FDR(%)
100
50
0 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34 0.36
Fig. 3. Tradeoffs between FAR and FDR of Monte Carlo runs from 20 to 50 s in two FD designs using uniform and Gaussian distributions for parametric uncertainties, respectively. 7. CONCLUSION For linear systems with polynomial dependence on probabilistic parametric uncertainties, this paper presents a new probabilistic robust parity relation approach using polynomial chaos theory. In contrast to existing polynomial chaos-based fault diagnosis literature, the proposed approach uses exact PCEs with no truncation errors, and is more computationally efficient, due to exploiting the polynomial parametric dependence within the system model.
Ackermann, J. (2002). Robust Control: The Parameter Space Approach. Springer-Verlag, London, UK. Blesa, J., Puig, V., and Saludes, J. (2012). Robust identification and fault diagnosis based on uncertain multiple input-multiple output linear parameter varying parity equations and zonotopes. Journal of Process Control, 22, 1890–1912. Blesa, J., Puig, V., Saludes, J., and Fern´andez-Cant´ı, R.M. (2016). Set-membership parity space approach for fault detection in linear uncertain dynamic systems. International Journal of Adaptive Control and Signal Processing, 30, 186–205. Budny, K. (2016). An extension of the multivariate Chebyshev’s inequality to a random vector with a singular covariance matrix. Communications in Statistics - Theory and Methods, 45(17), 5220–5223. Chen, X. and Zhou, K. (1997). On the probabilistic characterization of model uncertainty and robustness. In Proceedings of the 36th Conference on Decision and Control, 3616–3621. San Diego, USA. Chow, E.Y. and Willsky, A.S. (1984). Analytical redundancy and the design of robust failure detection systems. IEEE Transactions on Automatic Control, 29, 603–614. Ding, S.X. (2013). Model-Based Fault Diagnosis Techniques: Design Scheme, Algorithms, and Tools. Springer-Verlag, London, 2nd edition. Du, Y., Duever, T.A., and Budman, H. (2015). Fault detection and diagnosis with parametric uncertainty using generalized polynomial chaos. Computers and Chemical Engineering, 76, 63–75. Du, Y., Duever, T.A., and Budman, H. (2016). Generalized polynomial chaos based fault detection and classification for nonlinear dynamic processes. Industrial & Engineering Chemistry Research, 55, 2069–2082. Esfahani, P.M. and Lygeros, J. (2016). A tractable fault detection and isolation approach for nonlinear systems with probabilistic performance. IEEE Transactions on Automatic Control, 61, 633–647. Gertler, J. (1998). Fault Detection and Diagnosis in Engineering Systems. Dekker, New York. Maˆıtre, O.P.L. and Knio, O.M. (2010). Spectral Methods for Uncertainty Quantification: With Applications to Computational Fluid Dynamics. Springer Science & Business Media. Mesbah, A., Streif, S., Findeisen, R., and Braatz, R.D. (2014). Active fault diagnosis for nonlinear systems with probabilistic uncertainties. In Proceedings of the 19th IFAC World Congress, 7079–7084. Cape Town, South Africa. Ploix, S. and Adrot, O. (2006). Parity relations for linear uncertain dynamic systems. Automatica, 42, 1553–1562. Xiu, D. (2010). Numerical Methods for Stochastic Computations: A Spectral Method approach. Princeton University Press. Zhong, M., Zhang, L., Ding, S.X., and Zhou, D.H. (2016). A probabilistic approach to robust fault detection for a class of nonlinear systems. IEEE Transactions on Industrial Electronics. Doi: 10.1109/TIE.2016.2637308.
1047
ScienceDirect
IFAC PapersOnLine 50-1 (2017) 1019–1024 Probabilistic robust parity relation for fault Probabilistic robust parity relation for fault Probabilistic relation for fault detectionrobust using parity polynomial chaos detection using polynomial chaos detection using polynomial chaos ∗ ∗ ∗
Yiming Wan ∗ Eranda Harinath ∗ Richard D. Braatz ∗ Yiming Yiming Wan Wan ∗∗ Eranda Eranda Harinath Harinath ∗∗ Richard Richard D. D. Braatz Braatz ∗∗ ∗ Yiming Wan Eranda Harinath Richard D. Braatz Massachusetts Institute of Technology, Cambridge, MA, USA ∗ ∗ Massachusetts Institute of Technology, Cambridge, MA, USA ∗ Massachusetts Institute of Technology, Cambridge, MA, USA Massachusetts Institute Technology, Cambridge, MA, USA e-mail: {ywan,oferanda, braatz}@mit.edu e-mail: e-mail: {ywan, {ywan, eranda, eranda, braatz}@mit.edu braatz}@mit.edu e-mail: {ywan, eranda, braatz}@mit.edu Abstract: In this paper, a robust parity relation approach is proposed for fault detection Abstract: In paper, robust relation is for detection Abstract: In this thiswith paper, aprobabilistic robust parity parity relation approach approach is proposed proposed for fault fault detection of linear systems the a time-invariant parameter uncertainties. To deal with Abstract: In thiswith paper, aprobabilistic robust parity relation approach is proposed for fault detection of linear systems the time-invariant parameter uncertainties. To deal of linear systems with the time-invariant deal with with polynomial dependence on probabilistic uncertain parameters, a setparameter of parityuncertainties. relations withTopolynomial of linear systems with the probabilistic time-invariant parameter uncertainties. To deal with polynomial dependence on uncertain aa set of parity polynomial dependence on to uncertain parameters, set of condition. parity relations relations with polynomial parameterization is derived decoupleparameters, the unknown initial Due towith the polynomial polynomial dependence on to uncertain parameters, a set of condition. parity relations with polynomial parameterization is derived decouple the unknown initial Due to the parameterization derived to decouple the unknown Due to the polynomial polynomial structure having aisfinite degree, the generated residual initial vector condition. admits an exact polynomial chaos parameterization isfinite derived to decouple the unknown initial condition. Due to the polynomial structure having a degree, the generated residual vector admits an exact polynomial chaos structure a finite generated vectoronline admits an exact polynomial chaos expansionhaving with the samedegree, degree,the which enablesresidual an efficient computation for quantifying structure having a finite degree, the generated residual vector admits an exact polynomial chaos expansion with degree, which an efficient computation for expansion with the same degree, which enables an moments efficient online online computation foraquantifying quantifying the moments of the the same residual vector. Theenables obtained are used to calculate confidence expansion with the same degree, which enables an efficient online computation for quantifying the moments of the residual The obtained moments used calculate aa confidence the moments of the the consistency residual vector. vector. Thethe obtained moments are used to to calculate confidence set for checking between observed system are behavior and the parity relations. the moments of the the residual vector. Thethe obtained moments are used to calculate a confidence set for between observed system and relations. set for checking checking the consistency consistency between therelation observed system behavior behavior and the the parity parity relations. Compared to deterministic robust parity methods, the employment of probabilistic set for checkingdeterministic the consistency between therelation observed system behavior and the parity relations. Compared robust the employment of Compared toallows deterministic robust parity parity methods, the alarm employment of probabilistic probabilistic informationto a less conservative way relation to ensuremethods, a low false rate while maintaining Compared to deterministic robust parity relation methods, the employment of probabilistic information a conservative to low rate maintaining information allows a less less conservative way to ensure ensure low false false alarm alarm rate while while maintaining a high fault allows detection rate. In contrastway to the existingaa polynomial chaos-based fault diagnosis information allows a less conservative way the to ensure a polynomial low false alarm rate while maintaining aaliterature, high detection In chaos-based fault high fault fault detection rate. In contrast contrast to the existing existing polynomial chaos-based fault diagnosis diagnosis our proposedrate. approach avoidsto truncation errors of polynomial chaos expansions, and aliterature, high fault detection rate. In contrast to the existing polynomial chaos-basedexpansions, fault diagnosis our approach errors polynomial and literature, our proposed proposed demanding. approach avoids avoids truncation errors of of is polynomial chaos expansions, and is less computationally The truncation proposed approach illustratedchaos by using a two-tank literature, our proposed demanding. approach avoids truncation errors of polynomial chaosusing expansions, and is less is less computationally computationally demanding. The The proposed proposed approach approach is is illustrated illustrated by by using aa two-tank two-tank system example. is less computationally demanding. The proposed approach is illustrated by using a two-tank system system example. example. system © 2017, example. IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Fault detection, parametric uncertainties, parity relation, polynomial chaos. Keywords: Keywords: Fault Fault detection, detection, parametric parametric uncertainties, uncertainties, parity parity relation, relation, polynomial polynomial chaos. chaos. Keywords: Fault detection, parametric uncertainties, parity relation, polynomial chaos. 1. INTRODUCTION have a low probability of occurrence. For tractable online 1. have of tractable 1. INTRODUCTION INTRODUCTION have aa low low probability probability of occurrence. occurrence. For tractable online computation, a convex set needs to be For adopted as anonline outer 1. INTRODUCTION have a low probability of occurrence. For tractable online computation, a convex set needs to be adopted as computation, a to convex set needs adopted as an an outer outer The parity relation approach, initiated by Chow and approximation the exact set to of be system behaviors precomputation, a to convex exact set needs to be adopted as an outer The parity relation approach, by set of behaviors preThe parity relation approach,asinitiated initiated by Chow Chow and approximation to the the exact setrelations, of system system behaviors preWillsky (1984), is recognized one of main classesand of approximation dicted by the uncertain parity which introduces The parity relation approach,asinitiated by Chow and approximation to the exact setrelations, of system behaviors preWillsky (1984), is one of uncertain parity which Willsky (1984), is recognized recognized asdetection one of of main main classes of dicted dicted by by the the that uncertain parity relations, which introduces introduces techniques for model-based fault (FD)classes (Gertler, conservatism further reduces FD performance. Willsky (1984),model-based is recognized asdetection one of main classes of dicted by the uncertain parity relations, which introduces techniques (FD) techniques for model-based fault detection (FD) (Gertler, (Gertler, conservatism that that further further reduces reduces FD FD performance. performance. 1998; Ding,for 2013). The parityfault relation is a receding horizon conservatism techniques for model-based fault detection (FD) (Gertler, conservatism that further of reduces FD performance. By explicit consideration the probabilistic 1998; Ding, 2013). The parity relation is a receding horizon 1998; Ding, 2013). parity relation isfrom a receding horizon By explicit consideration of the probabilistic information input-output (I/O)The model decoupled the unknown information 1998; Ding, 2013). The parity relation is a receding horizon By uncertainties, explicit consideration of the probabilistic information the probabilistic robust FD approach input-output (I/O) decoupled the unknown input-output (I/O) model decoupled from the or unknown initial states by usingmodel a so-called parityfrom matrix vector. of By explicit consideration of the probabilistic information of uncertainties, the probabilistic robust FD approach input-output (I/O) model decoupled from the or unknown of uncertainties, the probabilistic robust FD approach has received some attention for reducing the conservatism initial states by using a so-called parity matrix vector. initial by using so-called the parity matrix orbetween vector. of Faults states are detected bya checking consistency uncertainties, the probabilistic robust FD approach has received some attention for reducing the conservatism initial states by using a so-called parity matrix or vector. has received some attention for reducing the conservatism of the deterministic robust approaches (Esfahani and Faults are detected by checking the consistency between Faults are relation detectedand bythe checking the system consistency between has the parity measured I/O behavior. received some attention for reducing the conservatism the deterministic robust approaches (Esfahani and Faults are relation detectedand bythe checking the system consistency between of of the deterministic robust approaches (Esfahani and Lygeros, 2016; Zhong et al., 2016). As a computational tool the parity measured I/O behavior. the parity relationfor and theinmeasured system I/O behavior. A key challenge FD real applications is to ensure Lygeros, of the deterministic robust approaches (Esfahani tool and 2016; Zhong et al., 2016). As a computational the parity relation and the measured system I/O behavior. Lygeros, 2016; Zhong et al., 2016). As a computational tool for probabilistic uncertainty propagation, the generalized A key for FD to A key challenge challenge for uncertainties FD in in real real applications applications isalarm to ensure ensure robustness to model for low falseis rate Lygeros, 2016; Zhong et al., 2016). As a computational tool for probabilistic uncertainty propagation, the generalized A key challenge for uncertainties FD in real applications isalarm to ensure for probabilistic the employed generalized chaosuncertainty (PC) theorypropagation, (Xiu, 2010) was in robustness tomaintaining model for low false falserate rate polynomial robustness model uncertainties low alarm rate (FAR) whileto high faultfor detection (FDR). probabilistic uncertainty propagation, the employed generalized polynomial chaos (PC) theory (Xiu, 2010) was in robustness tomaintaining model uncertainties for low falserate alarm rate for polynomial chaos (PC) theory (Xiu, 2010) was employed in the probabilistic robust FD methods proposed by Mesbah (FAR) while high fault detection (FDR). (FAR) whiledevelopment maintaining of high detection rateapproach (FDR). polynomial The early thefault parity relation chaosrobust (PC) theory (Xiu, 2010) was employed in the probabilistic FD methods proposed by Mesbah (FAR) while maintaining high fault detection rate (FDR). the probabilistic robust FD methods proposed by Mesbah et al. (2014); Du et al. (2015, 2016). By expressing the The early development of the parity relation approach The early development of the parity relation approach mainly focused on additive unknown disturbances, and the probabilistic robust FD methods proposed by Mesbah et al. (2014); Du et al. (2015, 2016). By expressing the The early development of the parity relation approach et al. (2014); Du et al. function (2015, 2016). By expressing the dependence as a series of orthogonal mainly focused on unknown and mainly focused on additive additive unknown disturbances, disturbances, and uncertainty could not effectively address multiplicative uncertainties. al. (2014); Du et al. function (2015, 2016). By expressing the uncertainty dependence as aa series of orthogonal mainly focused on additive unknown disturbances, and et uncertainty dependence function as series of orthogonal polynomials, the PC framework captures the evolution of could not effectively address multiplicative uncertainties. could not effectively address multiplicative uncertainties. uncertainty dependence function as a series of orthogonal polynomials, the PC PC framework framework captures the evolution evolution of could not effectively address multiplicative uncertainties. the captures the of By employing set-based tools, the set-membership parity polynomials, probability distributions with the PC expansion (PCE) polynomials, the PC framework captures the evolution of By employing set-based tools, parity probability These distributions with the PC expansion expansion (PCE) By employing set-based tools, the the set-membership parity coefficients. probability distributions with PC (PCE) relation approach was proposed toset-membership deal with determinisPC-based FDthe methods rely on applying By employing set-based tools, the set-membership parity probability distributions with the PC expansion (PCE) relation approach was proposed to deal with determiniscoefficients. These PC-based FD methods rely on applying relation approach was proposed to deal with determiniscoefficients. PC-based FD methods rely on applying tic unknown-but-bounded parametric uncertainties (Ploix the intrusiveThese Galerkin projection to the differential equarelation approach was proposed to deal with determiniscoefficients. These PC-based FD methods rely on applying tic uncertainties (Ploix Galerkin projection to the equatic unknown-but-bounded parametric uncertainties (Ploix the the intrusive intrusive Galerkin projection to which the differential differential equaandunknown-but-bounded Adrot, 2006; Blesa et parametric al., 2012, 2016). The approach tions describing the plant dynamics, limits their aptic unknown-but-bounded parametric uncertainties (Ploix the intrusive Galerkin projection to which the differential equaand Adrot, 2006; Blesa et al., 2012, 2016). The approach tions describing the plant dynamics, limits their apand Adrot, 2006; Blesa et al., 2012, 2016). The approach tions describing the plant dynamics, which limits their apdecouples the unknown initial states by introducing the plicability to stable systems only, because the convergence and Adrot,the 2006; Blesa et al., 2012, 2016). The approach tions describing the systems plant dynamics, whichthe limits their apdecouples unknown initial by the stable because decouples the as unknown initial states by introducing introducing the plicability plicability todoes stable only, because the convergence convergence parity matrix the function of states the uncertain parameters. of the PCEto notsystems hold foronly, unstable systems. Moreover, decouples the as unknown initial states by introducing the plicability todoes stable systems only, because the convergence parity matrix function of uncertain parameters. PCE not for Moreover, parity matrix in as the the function of the the uncertain parameters. of of the the PCE does to notbehold hold for unstable unstable systems. Moreover, This results a set of parity relations parameterized the PCE needs truncated, and systems. no systematic way parity matrix in as the function of the uncertain parameters. of the PCE does to notbehold for unstable systems. Moreover, This results a set of parity relations parameterized the PCE needs truncated, and no systematic This in aparameters, set of parity relations theavailable PCE needs to be truncated, and no systematic way by theresults uncertain rather than aparameterized single parity is to account for such truncation errors in way the This results in aparameters, set of parity relations parameterized theavailable PCE needs to be truncated, and no systematic way by the rather than parity to for by the uncertain uncertain rather than aa single single parity is available to account account for such such truncation truncation errors errors in in the the relation as in the parameters, case of additive disturbances. Then the is robust FD system design. by the uncertain parameters, rather than a single parity availablesystem to account for such truncation errors in the relation as case of Then the is robust relation as in in isthe the case by of additive additive disturbances. the robust FD FD system design. design. FD decision made checkingdisturbances. whether the Then observed relation as in isthe case by of additive disturbances. Then the robustwork FD system design.the robust FD problem for lininvestigates FD decision whether the observed FD decision is made made by checking checking whether observed system behavior is admitted by the set ofthe parity rela- This Thissystems work investigates investigates the robust robust FD problem problem for linlinFD decision is made by checking whether the observed This work the FD for ear with polynomial dependence on probabilissystem behavior is admitted by the set of parity relasystem behavior such is admitted by the robust set of parity rela- This tions. Although a deterministic work investigates the robust FD problem for linear systems with polynomial dependence on probabilissystem behavior is admitted by the set of parity relaear systems with polynomial dependence on probabilistic time-invariant parametric uncertainties. The above tions. Although such aa deterministic robust parity tions. Althoughcan such deterministic robust even parityinrelarelation approach ensure zero false alarms the ear systems with polynomial dependence on The probabilistic parametric above tions. Although such a deterministic robust even parity relatic time-invariant time-invariant parametric uncertainties. uncertainties. The robust above of the deterministic and probabilistic tion ensure false tion approach can ensure zero false alarms alarms even in in the the limitations worstapproach case, the can resulting FDzero performance is compromised, tic time-invariant parametric uncertainties. The above limitations of the deterministic and probabilistic robust tion approach can ensure zero false alarms even in the FD limitations of the deterministic and probabilistic robust methods can be effectively handled by our proposed worst case, the resulting FD performance is compromised, worst case, the worst-case resulting FDparameter performance is compromised, because some uncertainties often limitations of the deterministic and probabilistic robust can effectively handled by proposed worst case, the worst-case resulting FDparameter performance is compromised, FD methods methods robust can be be parity effectively handled by our ourusing proposed probabilistic relation approach PC. because some uncertainties often FD because some worst-case parameter uncertainties often FD methods robust can be parity effectively handled by ourusing proposed ⋆ probabilistic relation approach PC. Support some from Novartis is acknowledged. because worst-case parameter uncertainties often probabilistic robust parity relation approach using PC. ⋆ probabilistic robust parity relation approach using PC. ⋆ Support Support from from Novartis Novartis is is acknowledged. acknowledged.
⋆ Support from Novartis is acknowledged. Copyright © 2017 IFAC 1042 2405-8963 © 2017, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Copyright © 2017 IFAC 1042 Copyright ©under 2017 responsibility IFAC 1042Control. Peer review of International Federation of Automatic Copyright © 2017 IFAC 1042 10.1016/j.ifacol.2017.08.210
Proceedings of the 20th IFAC World Congress 1020 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
Similarly to Ploix and Adrot (2006); Blesa et al. (2016), a parity matrix is designed as a polynomial function of the uncertain parameters. By exploiting the polynomial parametric dependence, both the uncertain parity relation and the generated residual vector are written into a linear form with respect to the monomials of uncertain parameters, which enables efficient online stochastic uncertainty quantification using offline PCE of the monomials. Moreover, due to the polynomial parametric dependence, the involved PCEs are exact without any truncation errors. Using the moment information computed from the PCE of the residual vector, a Chebyshev bound is calculated for checking the consistency between the parity relation and the observed I/O data. This paper is organized as follows. Section 2 states the probabilistic FD problem. Section 3 presents the computation of the polynomial parity matrix. The residual uncertainty quantification is discussed in Section 4. Section 5 constructs the consistency test for the FD decision making. The simulation results are presented in Section 6. 2. PROBLEM STATEMENT 2.1 System description Consider the linear uncertain dynamical system x(k + 1) = A(θ)x(k) + B(θ)u(k) y(k) = C(θ)x(k) + D(θ)u(k) + e(k)
(1)
where the system inputs u(k) ∈ Rnu and outputs y(k) ∈ Rny are measurable, the states x(k) ∈ Rnx are unmeasured, e(k) denotes the additive stochastic measurement noises with zero mean and known covariance Σe , and A(θ), B(θ), C(θ), and D(θ) are matrices whose entries are polynomials of the uncertain parameters θ ∈ Rnθ . Such parametric uncertainties are often encountered in practice when specific knowledge is available on how uncertain parameters enter a system model (Ackermann, 2002). The uncertain parameters θ are time invariant, which is also a good approximation for slowly time-varying parameters (Blesa et al., 2016). We also assume that the entries of θ are mutually independent random variables with known time-invariant probability density functions (PDFs) and supports. Such a description of θ can be derived from the a prior knowledge of the underlying system, or offline parameter identification from data.
yk,h
y(k − h + 1) y(k − h + 2) , = .. .
C(θ) C(θ)A(θ) Oh (θ) = .. .
,
C(θ)Ah−1 (θ) 0 ··· 0 . .. . .. C(θ)B(θ) D(θ) Hh (θ) = .. .. .. . 0 . . C(θ)Ah−2 (θ)B(θ) C(θ)Ah−3 (θ)B(θ) · · · D(θ) (3) and uk,h and ek,h are defined similarly to yk,h . A parity matrix V (θ) is designed to satisfy V (θ)Oh (θ) = 0 for any θ (4) so that the premultiplication of (2) by this matrix generates the parity relation (5) V (θ)(yk,h − Hh (θ)uk,h − ek,h ) = 0, without the unknown initial term Oh (θ)x(k − h + 1). The parity matrix V (θ) is assumed to be a polynomial matrix function of the uncertain parameters θ, due to the polynomial dependence of Oh (θ) on θ. Using this parity relation (5), the residual vector is generated with polynomial dependence on θ and e(k): r(k, θ, e) = G(θ)zk,h − V (θ)ek,h (6) [ ⊤ ] ⊤ ⊤ with G(θ) = [V (θ) −V (θ)Hh (θ)] and zk,h = yk,h uk,h .
y(k) D(θ)
With the bounded descriptions θ ∈ Θ and e(k) ∈ E, the set-based consistency test was proposed as (Ploix and Adrot, 2006; Blesa et al., 2016) (7) 0 ∈ Γk , where Γk = {r|r = G(θ)zk,h − V (θ)ek,h , θ ∈ Θ, e(k) ∈ E}. When (7) holds, there exists θ0 ∈ Θ and {e0 (i) ∈ E}ki=k−h+1 such that r(k, θ0 , e0 ) = 0 and the occurrence of faults cannot be concluded. Otherwise, the observed I/O data are inconsistent with the parity relation, which indicates the presence of faults. The exact computation of the set Γk in (7) is often prohibitive due to its nonlinear dependence on θ. This difficulty was circumvented in Blesa et al. (2016) by adopting a different linear parametrization, at the cost of reduced FDR. Since the residual r(k, θ, e) is a polynomial vector of θ and e, (6) was equivalently rewritten as (Blesa et al., 2016) r (k, ξ(θ), n(θ, ek,h )) = η(k, ξ(θ)) + n(θ, ek,h ), (8) ⊤
where ξ(θ) = [ξ0 (θ) · · · ξNG −1 (θ)] , η(k, ξ(θ)) =
2.2 Robust parity relation
N∑ G −1
φi (zk,h )ξi (θ),
(9)
i=0
With the unknown-but-bounded assumption on the parametric uncertainties θ and the measurement noises e(k) in the system (1), Ploix and Adrot (2006) proposed a deterministic robust parity relation approach which is briefly reviewed below to provide background and to motivate our investigated problem. From the system model (1), the stacked output equation over a time window [k − h + 1, k] can be derived as where
yk,h = Oh (θ)x(k − h + 1) + Hh (θ)uk,h + ek,h
(2)
(10) n(θ, ek,h ) = −V (θ)ek,h , NG represents the number of monomials in the polynomial matrix G(θ), ξi (θ) denotes the monomial basis, and the corresponding monomial coefficient φi (zk,h ) is linear in the I/O data zk,h . With the new parameters ξ and n, the set Γk in (7) is replaced by (11) Υk = {r | r = η(k, ξ) + n, ξ ∈ Ξ, n ∈ N } . For the sake of tractable computation, the sets Ξ and N are used as outer approximations of {ξ|ξ = ξ(θ), θ ∈ Θ} and {n|n = V (θ)ek,h , θ ∈ Θ, e(k) ∈ E}. Therefore, Γk ⊆ Υk , and using 0 ∈ Υk as the consistency test reduces FDR.
1043
Proceedings of the 20th IFAC World Congress Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
Being a worst-case formulation, the above deterministic parity relation approach does not take into account the probability of occurrence for all possible values of parametric uncertainties and noises, even for values that have vanishingly low probability (e.g., could occur once in 10100 years). As a result, the deterministic robust approach is conservative in ensuring FAR, thus compromising the FDR. This article employs PC theory in the robust parity relation so that probability information of parametric uncertainties in the linear system (1) can be incorporated to enhance FDR while ensuring a low FAR. 3. COMPUTING A POLYNOMIAL PARITY MATRIX This section presents how to compute the polynomial parity matrix V (θ). For the ease of understanding, we discuss first the case of a univariate uncertain parameter θ ∈ R , and then proceed to the general case of the multivariate uncertain parameters θ ∈ Rnθ . For a scalar θ, assume that the degrees of the polynomial matrices A(θ) and C(θ) are dA and dC , respectively. Then the observability matrix Oh (θ) has a degree of dO = dC + (h − 1)dA . (12) When a polynomial parity matrix V (θ) is selected with a degree dV , (4) becomes ) d (d V O ∑ ∑ Vi θ i Oj θj = 0. (13) i=0
j=0
The multiplication on the left-hand side of (13) is a polynomial with a degree dV + dO + 1, whose coefficient ∑i of the monomial basis with a degree i is j=0 Vj Oi−j . All of these monomial coefficients need to be zero, i.e., ∑i j=0 Vj Oi−j = 0, so that (13) holds for all θ. Equation 13 can be compactly written as O 0 O1 O2 · · · 0 0 O0 O1 · · · 0 0 0 O ··· 0 = 0. 0 (14) [V0 V1 · · · VdV ] �� � � ... ... ... . . . ... V 0 0 0 · · · OdO � �� � M
The multiplication of V and the ith block-column of M is the coefficient of the monomial basis with a degree i − 1. The maximum number of block-columns in M is determined by the degree of V (θ)Oh (θ), and the size of M is (dV + 1)hny × nM , with nM ≤ (dV + dO + 1)nx . This idea can be extended to the case of a nθ -dimensional ⊤ parameter vector θ = [θ1 θ2 · · · θnθ ] . Let θ[S] represent ∑nθ snθ s1 s2 a monomial basis θ1 θ2 · · · θnθ with a degree i=1 si and S = (s1 , s2 , · · · , snθ ). Then a ∑ polynomial matrix M (θ) can be decomposed as M (θ) = M[Si ] θ[Si ] , where M[Si ] is Si ∈S
the coefficient matrix for the monomial basis θ[Si ] , and the set S is determined by all the monomials presented in the polynomial matrix M (θ). If the degree of M (θ) is dM , the maximum number of monomial terms in M (θ) is (dM + nθ )! . (15) g(dM , nθ ) = dM !nθ !
1021
Based on the definitions, (4) becomes ( ) ∑ ∑ V[Si ] θ[Si ] O[Qi ] θ[Qi ] = 0 Si ∈S
(16)
Qi ∈Q
for multivariate uncertain parameters. Finding a polynomial parity matrix V (θ) again relies on solving VM = 0, with V and M defined similarly to (14). Now the size of M is g(dV , nθ )hny × nM , with nM ≤ g(dV + dO , nθ )nx . The number of block-columns nM can be smaller than g(dV +dO , nθ )nx due to the possible zero coefficients O[Qi ] . The existence of a solution to VM = 0 is ensured by having the matrix M have more rows than columns, i.e., g(dV , nθ )hny > g(dV + dO , nθ )nx . Using (15), this inequality can be rewritten as ) nθ ( ∑ dO nx > 1− , (17) d + d + i hn V O y i=1
which implies that by selecting a sufficiently high degree dV and a large length h of the time window, a polynomial parity matrix V (θ) exists so that (14) holds. 4. QUANTIFYING RESIDUAL UNCERTAINTY With the parity matrix obtained in Section 3, the residual vector can be generated from (6). As a function of the random parameters θ and the stochastic noises e(k), the generated residual vector is also stochastic. In this section, PC theory is applied for efficient online uncertainty quantification of the residual vector. 4.1 PCE of multivariate polynomials Consider a multivariate polynomial matrix denoted by m ∑ M (θ) = Mi ξi (θ), (18) i=0
with ξi (θ) representing a monomial basis with its degree at most p. Then M (θ) has an exact PCE representation M (θ) =
Np ∑
Ψj Φj (θ),
(19)
j=0
where Φj (θ) denotes the multivariate PC basis function. The total number of terms in (19) is Np + 1 = (nnθθ+p)! !p! , depending on the dimension nθ of θ and the order p of M (θ) (Xiu, 2010). The PCE coefficients {Ψj } for M (θ) can be efficiently computed by following two steps: (i) determine the PCE coefficients {βij } of each monomial basis ξi (θ), i.e., ξi (θ) =
Np ∑
βij Φj (θ);
(20)
j=0
and (ii) compute the PCE coefficients {Ψj } in (19) as m ∑ Mi βij (21) Ψj = i=0
by substituting (20) into (18). Example 1. The Legendre polynomials form an orthogonal basis with respect to a uniform random scalar θ over
1044
Proceedings of the 20th IFAC World Congress 1022 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
[−1, 1]. The first four Legendre polynomials are given by (Maˆıtre and Knio, 2010) 1 1 0 0 0 Φ0 (θ) Φ1 (θ) 0 1 0 0 θ Φ (θ) = − 1 0 3 0 θ2 . 2 2 2 Φ3 (θ) 0 − 32 0 52 θ3
By inverting the coefficient matrix, the monomials θ, θ2 , and θ3 can be represented as Φ0 (θ) θ 0 1 0 0 Φ1 (θ) θ 2 = 1 0 2 0 Φ (θ) , 3 3 2 3 3 2 0 5 0 5 θ Φ (θ) 3
where the elements of each row of the coefficient matrix are the PCE coefficients of the corresponding monomial. Then the PCE coefficients of any polynomial matrix (18) with a degree less than 4 can be determined according to (21).
The PC basis functions belong to the Askey scheme of orthogonal polynomials, satisfying (Xiu, 2010) ⟨Φi (θ), Φj (θ)⟩ = ⟨Φ2i (θ)⟩δij , where δij is the Kronecker delta, and ⟨·, ·⟩ defines the inner product ∫ ⟨α(θ), β(θ)⟩ = α(θ)β(θ)fθ dθ
respectively. The polynomial dependence of η(k, ξ(θ)) and n(θ, ek,h ) on θ in (24) allows the uncertainty quantification to be decomposed into four sequential steps: (i) For the selected polynomial bases {Φj (θ)}, compute the inner products {⟨Φ2j (θ)⟩} and the PCE coefficients of each monomial basis ξi (θ) as in (20). (ii) Compute Cov{n(θ, ek,h )}. From (10), it follows that ⊤ Cov{n(θ, ek,h )} = E{V (θ)ek,h e⊤ k,h V (θ)}
(25) = E{V (θ)Σe,h V ⊤ (θ)} = E{Σn (θ)}, where Σe,h = Ih ⊗ Σe denotes the covariance of ek,h . Note that Σn (θ) = V (θ)Σe,h V ⊤ (θ) in (25) has polynomial dependence on θ. After calculating the PCE coefficients of Σn (θ) from (21), E{Σn (θ)} is determined by its zero-order PCE coefficient according to (22). (iii) Compute E{η(k, ξ(θ))} and Cov{η(k, ξ(θ))} by using the online data zk,h . For a given time window, the I/O data zk,h can be regarded as deterministic. Then, by substituting (20) into (9), the PCE of η(k, ξ(θ)) is η(k, ξ(θ)) =
E{M (θ)} =
Np ∑ j=0
Cov{M (θ)} =
Ψj ⟨Φj (θ), 1⟩ = Ψ0 ,
N p Np ∑ ∑ i=1 j=1
=
Np ∑ j=1
(22)
Ψi Ψ⊤ j ⟨Φi (θ), Φj (θ)⟩
2 Ψj Ψ⊤ j ⟨Φj (θ)⟩.
(23)
Note that the inner products of basis functions {⟨Φ2j (θ)⟩} 1 for an uniform can be pre-computed, e.g., ⟨Φ2j (θ)⟩ = 2j+1 random scalar θ, so that the covariance in (23) can be efficiently calculated using PCE coefficients. 4.2 Polynomial chaos approach to computing residual mean and covariance In this subsection, the mean and covariance of the residual vector (6) is computed by applying PC theory reviewed in Section 4.1. For this purpose, the linear parameterization form (8), instead of (6), is adopted. Due to the probabilistic independence between θ and e(k), E{n(θ, ek,h )} = 0 and E{n(θ, ek,h )η ⊤ (k, ξ(θ))} = 0 are obtained, thus the mean and covariance of the uncertain residual vector (8) are µk = E{η(k, ξ(θ))}, and (24) Σk = Cov{η(k, ξ(θ))} + Cov{n(θ, ek,h )},
ηj (k)Φj (θ),
j=0
Ω
with respect to the PDF fθ of θ over the support Ω. The selection of PC bases is determined by probability distributions of random parameters, e.g., Legendre polynomials for uniform distributions, Hermite polynomials for Gaussian distributions. Exploiting the above orthogonality, the mean and covariance of the polynomial matrix M (θ) can be conveniently determined from its PCE coefficients {Ψj } by (Xiu, 2010)
Nq ∑
ηj (k) =
N∑ G −1
(26) φi (zk,h )βij ,
i=0
with q representing the degree of η(k, ξ(θ)). The PCE coefficients ηj (k) can be efficiently computed by using the online data zk,h and the PC coefficients βij obtained offline in (20). From the PCE coefficients ηj (k) in (26), E{η(k, ξ(θ))} and Cov{η(k, ξ(θ))} can be calculated from (22) and (23). (iv) According to (24), compute the mean µk and covariance Σk of the residual vector (8) using results from Steps (ii) and (iii). The first two steps can be performed offline, and the last two steps are online computations involving only basic matrix manipulations, which results in much more efficient computation, compared to using Galerkin projection based methods in Mesbah et al. (2014); Du et al. (2015, 2016). The involved PCEs are all exact with no truncation errors, due to the polynomial parametric dependence of a finite order in the system model (1). In contrast, the PCE-based methods in Mesbah et al. (2014); Du et al. (2015, 2016) have truncation errors that are not explicitly considered. 5. CONSISTENCY TEST In this section, the Chebyshev bound of the residual vector are calculated for the consistency test by using the residual mean µk and covariance Σk obtained in Section 4.2. Lemma 1. (Chen and Zhou (1997); Budny (2016)). For the stochastic residual vector rk = r(k, θ, e) with the mean µk and the covariance Σk ,
1045
(i) Assume Σk is positive definite with rank(Σk ) = nr . For any ϵ > 0, it holds that { } nr ⊤ , Pr (rk − µk ) Σ−1 k (rk − µk ) < ϵ > 1 − ϵ
Proceedings of the 20th IFAC World Congress Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
(ii) Assume Σk is rank deficient with rank(Σk ) = mk , and has a eigenvalue decomposition Σk = Uk Πk Uk⊤ , where the diagonal matrix Πk = diag(λ1 , · · · , λmk , 0, · · · , 0), with λi > 0. For any ϵ > 0, it holds that } { mk ⊤ , Pr (rk − µk ) Uk Tk Uk⊤ (rk − µk ) < ϵ > 1 − ϵ −1 where Tk = diag(λ−1 1 , · · · , λmk , 0, · · · , 0). For a nonsingular or singular covariance matrix Σk , Lemma 1 provides a confidence ellipsoid ∆k with a confidence level 1 − γ, i.e., Pr {rk ∈ ∆k } > 1 − γ, (27) where ∆k and γ are defined by { � } nr � ⊤ , or ∆k = rk �(rk − µk ) Σ−1 k (rk − µk ) < ϵ , γ = { � } ϵ m � k ⊤ , ∆k = rk �(rk − µk ) Uk Tk Uk⊤ (rk − µk ) < ϵ , γ = ϵ respectively. With a large ϵ, rk ∈ ∆k holds with a high confidence level 1 − γ, then the ellipsoid ∆k describes the most plausible region in which the residual vector lies. For this reason, the confidence ellipsoid ∆k can be used in the consistency test as Γk in (7), i.e., performing the consistency test by checking whether 0 ∈ ∆k holds. This approach is as simple as checking whether J < 1, with ⊤
J = ϵ−1 (rk − µk ) Σ−1 k (rk − µk ) , or
(28) ⊤ J = ϵ−1 (rk − µk ) Uk Tk Uk⊤ (rk − µk ) . The invalidation of J < 1 indicates the presence of faults. Remark 1. Since Lemma 1 holds for any PDF with the given mean and covariance, the confidence level 1 − γ in (27) is always conservative, i.e., the actual probability of rk ∈ ∆k can be much larger than 1 − γ. Therefore, it is possible to increase γ for an improved FDR without reducing FAR. 6. SIMULATION EXAMPLE Consider a two-tank process described by (Blesa et al., 2016) dx1 α1 √ α2 √ κ1 =− 2gx1 + 2gx2 + v1 dt A1 A1 A1 (29) dx2 α2 √ κ2 =− 2gx2 + v2 dt A2 A2 where α1 = α2 = 7 cm2 are the cross-sections of the two outlet holes, A1 = A2 = 28 cm2 are the cross-sections of the two tanks, κ1 = 2.33 cm3 /V s and κ2 = 1.34 cm3 /V s are the coefficients of the two inlet flows, v1 and v2 are input voltages to the pumps, x1 and x2 are the water levels of the two tanks, and g = 981 cm/s2 . The sampled ⊤ system inputs and outputs are u(k) = [v1 (k) v2 (k)] and ⊤ y(k) = [x1 (k) x2 (k)] + e(k), with the sampling interval ts = 1 s and the zero-mean white noises whose covariance is Σe = diag(10−4 , 10−4 ) cm2 . The water levels may slowly vary within the interval [xci − ∆xi , xci + ∆xi ], with xc1 = 30 cm, xc2 = 15 cm, and ∆x1 = ∆x2 = 5 cm. The continuous-time system (29) is discretized by the firstorder Euler method and linearized around the equilibrium point described by x01 , x02 , v10 , v20 , which results in the discrete-time model,
[
1023
] [ ] ][ δx1 (k + 1) a (x0 (k)) a12 (x02 (k)) δx1 (k) = 11 1 δx2 (k + 1) 0 a22 (x02 (k)) δx2 (k) ][ ] [ (30) δv1 (k) b11 0 , + 0 b22 δv2 (k)
where δvi = vi −vi0 , δxi = xi −x0i , and aij and bij are determined by constants of the system (29) and the uncertain water levels. To apply the proposed method, the nonlinear functions a11 (x01 (k)), a12 (x02 (k)), and a22 (x02 (k)) over the water level variation intervals [xci − ∆xi , xci + ∆xi ] were approximated by their Taylor expansions at the centers xc1 and xc2 , i.e., a ˆ11 (∆x1 (k)), a ˆ12 (∆x2 (k)), and a ˆ22 (∆x2 (k)), respectively. For this system, second-order Taylor expansions are sufficiently accurate. By doing so, the model (30) has polynomial dependence on two uncertain parameters θ1 = ∆x1 and θ2 = ∆x2 . Due to lack of in-depth probabilistic information, θ1 and θ2 are both assumed to be uniformly distributed in [−5, 5] cm, thus Legendre polynomial bases are adopted. The time window length h of the parity relation is selected to be h = 2, and a 2-by-4 polynomial parity matrix with a degree 2 is found. For both fault-free and faulty scenarios, 100 Monte Carlo (MC) runs are implemented with different time-invariant parameters θ in each run, with the fault being a leak in tank 1 injected at 20 s by adding ∆α1 = 0.01 cm2 to α1 in (29). Confidence ellipsoids ∆k and residual vectors of 100 MC runs with this fault are visualized in Fig. 1 for a sequence of 4 sampling instants, with γ = 0.15. The consistency condition 0 ∈ ∆k becomes invalid after 21 s, which is prompt detection of the leakage fault. The ensemble of the residual evaluations (28) in the fault-free and fault MC runs for the same γ shows zero FAR and high FDR accounting for the variations of uncertain parameters and measurement noises (Fig. 2). With a smaller leakage fault ∆α1 = 0.008 cm2 , the case of θ1 and θ2 being Gaussian distributed is also simulated. Each uncertain parameter has its standard deviations equal to 1 cm, so that its value almost surely belongs to the 5-sigma interval [−5, 5] cm. Two FD designs are considered: (i) uniform distribution over [−5, 5] cm is used for θ1 and θ2 in the design; and (ii) zero-mean Gaussian distribution with its standard deviation being equal to 1 cm is used in the design. For statistical performance evaluation, the FAR and FDR over a finite time window [k1 , k2 ] was used as defined by ∑k 2
number of false alarms at time instant i (number of fault-free MC runs) × (k2 − k1 + 1) ∑k2 i=k1 number of FDs at time instant i . FDR = (number of fault MC runs) × (k2 − k1 + 1) FAR =
i=k1
As shown in Figure 3, increasing γ from 0.1 to 0.26 leads to an improved FDR without sacrificing FAR, due to the reason explained in Remark 1. Moreover, the design (ii) has much higher FDR than the design (i) for all values of γ ∈ [0.1, 0.36], while it only suffers from higher FAR for γ = 0.3. This improvement of FD performance in the design (ii) benefits from exploiting the more informative probabilistic distribution information on the uncertain parameters.
1046
Proceedings of the 20th IFAC World Congress 1024 Yiming Wan et al. / IFAC PapersOnLine 50-1 (2017) 1019–1024 Toulouse, France, July 9-14, 2017
t = 20s residuals
0
-0.05 -0.1
confidence ellipsoids
0
r(2) k
r(2) k
0.05
REFERENCES
-0.05 -0.05
0
0.05
-0.05
0
r(1) k
t = 23s
0
r(2) k
r(2) k
0.05
r(1) k
t = 22s
0
Future work will focus on comparisons with existing setmembership parity relation methods, and extensions to nonlinear dynamical systems.
t = 21s 0.05
-0.05
-0.05 -0.1
-0.1 -0.05
0
0.05
-0.05
0
r(1) k
0.05
r(1) k
1
Normalized evaluation 0.5
Consistency evaluation
Consistency evaluation
Fig. 1. Confidence ellipsoids and residual vectors of 100 fault Monte Carlo runs from 20 to 23 s, with γ = 0.15.
0 10
15
6
Normalized threshold
20
25
30
35
40
45
50
45
50
Fault Injection
4 2 0 10
15
20
25
30
Time (s)
35
40
Fig. 2. The ensemble of the residual evaluations (28) in the fault-free and fault Monte Carlo runs for γ = 0.15.
FAR(%)
15 10
Uniform distribution Gaussian distribution
5 0 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34 0.36
FDR(%)
100
50
0 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34 0.36
Fig. 3. Tradeoffs between FAR and FDR of Monte Carlo runs from 20 to 50 s in two FD designs using uniform and Gaussian distributions for parametric uncertainties, respectively. 7. CONCLUSION For linear systems with polynomial dependence on probabilistic parametric uncertainties, this paper presents a new probabilistic robust parity relation approach using polynomial chaos theory. In contrast to existing polynomial chaos-based fault diagnosis literature, the proposed approach uses exact PCEs with no truncation errors, and is more computationally efficient, due to exploiting the polynomial parametric dependence within the system model.
Ackermann, J. (2002). Robust Control: The Parameter Space Approach. Springer-Verlag, London, UK. Blesa, J., Puig, V., and Saludes, J. (2012). Robust identification and fault diagnosis based on uncertain multiple input-multiple output linear parameter varying parity equations and zonotopes. Journal of Process Control, 22, 1890–1912. Blesa, J., Puig, V., Saludes, J., and Fern´andez-Cant´ı, R.M. (2016). Set-membership parity space approach for fault detection in linear uncertain dynamic systems. International Journal of Adaptive Control and Signal Processing, 30, 186–205. Budny, K. (2016). An extension of the multivariate Chebyshev’s inequality to a random vector with a singular covariance matrix. Communications in Statistics - Theory and Methods, 45(17), 5220–5223. Chen, X. and Zhou, K. (1997). On the probabilistic characterization of model uncertainty and robustness. In Proceedings of the 36th Conference on Decision and Control, 3616–3621. San Diego, USA. Chow, E.Y. and Willsky, A.S. (1984). Analytical redundancy and the design of robust failure detection systems. IEEE Transactions on Automatic Control, 29, 603–614. Ding, S.X. (2013). Model-Based Fault Diagnosis Techniques: Design Scheme, Algorithms, and Tools. Springer-Verlag, London, 2nd edition. Du, Y., Duever, T.A., and Budman, H. (2015). Fault detection and diagnosis with parametric uncertainty using generalized polynomial chaos. Computers and Chemical Engineering, 76, 63–75. Du, Y., Duever, T.A., and Budman, H. (2016). Generalized polynomial chaos based fault detection and classification for nonlinear dynamic processes. Industrial & Engineering Chemistry Research, 55, 2069–2082. Esfahani, P.M. and Lygeros, J. (2016). A tractable fault detection and isolation approach for nonlinear systems with probabilistic performance. IEEE Transactions on Automatic Control, 61, 633–647. Gertler, J. (1998). Fault Detection and Diagnosis in Engineering Systems. Dekker, New York. Maˆıtre, O.P.L. and Knio, O.M. (2010). Spectral Methods for Uncertainty Quantification: With Applications to Computational Fluid Dynamics. Springer Science & Business Media. Mesbah, A., Streif, S., Findeisen, R., and Braatz, R.D. (2014). Active fault diagnosis for nonlinear systems with probabilistic uncertainties. In Proceedings of the 19th IFAC World Congress, 7079–7084. Cape Town, South Africa. Ploix, S. and Adrot, O. (2006). Parity relations for linear uncertain dynamic systems. Automatica, 42, 1553–1562. Xiu, D. (2010). Numerical Methods for Stochastic Computations: A Spectral Method approach. Princeton University Press. Zhong, M., Zhang, L., Ding, S.X., and Zhou, D.H. (2016). A probabilistic approach to robust fault detection for a class of nonlinear systems. IEEE Transactions on Industrial Electronics. Doi: 10.1109/TIE.2016.2637308.
1047