1

Sensor Fault Detection, Isolation and Identification
Using Multiple Model-based Hybrid Kalman Filter
for Gas Turbine Engines*

arXiv:1505.02063v2 [] 27 Aug 2015

Bahareh Pourbabaee1 , Nader Meskin2 and Khashayar Khorasani1

Abstract
In this paper, a novel sensor fault detection, isolation and identification (FDII) strategy is proposed
by using the multiple model (MM) approach. The scheme is based on multiple hybrid Kalman filters (HKF) which represents an integration of a nonlinear mathematical model of the system with a
number of piecewise linear (PWL) models. The proposed fault detection and isolation (FDI) scheme is
capable of detecting and isolating sensor faults during the entire operational regime of the system by
interpolating the PWL models using a Bayesian approach. Moreover, the proposed multiple HKF-based
FDI scheme is extended to identify the magnitude of a sensor fault by using a modified generalized
likelihood ratio (GLR) method which relies on the healthy operational mode of the system. To illustrate
the capabilities of our proposed FDII methodology, extensive simulation studies are conducted for a
nonlinear gas turbine engine. Various single and concurrent sensor fault scenarios are considered to
demonstrate the effectiveness of our proposed on-line hierarchical multiple HKF-based FDII scheme
under different flight modes. Finally, our proposed HKF-based FDI approach is compared with various
filtering methods such as the linear, extended, unscented and cubature Kalman filters (LKF, EKF, UKF
and CKF, respectively) corresponding to both interacting and non-interacting multiple model (MM)
based schemes. Our comparative studies confirm the superiority of our proposed HKF method in terms
of promptness of the fault detection, lower false alarm rates, as well as robustness with respect to the
engine health parameters degradations.

I. INTRODUCTION
An effective fault detection, isolation and identification (FDII) technology can play a crucial
role in improving the system availability, safety and reliability as well as reducing the maintenance costs and risks of catastrophic failures. Over the past few years, many researchers have
focused on proposing sophisticated fault detection and isolation (FDI) schemes constituting as
a significant component of an FDII solution [1] and [2]. The general methodology for an FDII
consists of several important steps including: (a) generation of residuals as indicators of faults,
(b) isolation of the faulty actuator, or sensor, or component element, and (c) identification and
estimation of the severity of the fault parameter.
The FDII approaches can be categorized into three distinct groups including the modelbased, data-based and hybrid methods [3] and [4]. The model-based FDI approaches have been
introduced in survey papers [1], [2] and [5]–[7] and can be grouped into basic approaches,
namely (i) Kalman filter and unknown input observer based approaches, (ii) parity relations,
(iii) optimization-based algorithms, and (iv) parameter estimation and identification techniques.
*This publication was made possible by NPRP grant No. 4 - 195 - 2 - 065 from the Qatar National Research Fund (a member
of Qatar Foundation). The statements made herein are solely the responsibility of the authors.
1
B. Pourbabaee and K. Khorasani are with the Department of Electrical and Computer Engineering, Concordia University,
Montreal, Canada. b pourba@ece.concordia.ca and kash@ece.concordia.ca.
2
N. Meskin is with the Department of Electrical Engineering, Qatar University, Doha, Qatar. nader.meskin@qu.edu.qa.

August 28, 2015

DRAFT

2

Multiple model approach is a popular structure for FDII. It enables one not only to detect and
isolate different faults but it also provides one with information on the magnitude of the fault
and its identification information [8]–[13].
Since early research on FDII, gas turbines have been one of the challenging application areas
which have received much attention. The main approach in gas turbines FDII is based on the
Gas Path Analysis (GPA), which enables one to accomplish the actuator, sensor and component
diagnosis by observing the engine’s parameters such as the rotor speed, temperature and pressure
at different stages and the fuel flow rates [14]. Various fault diagnosis techniques have been
developed for gas turbines based on GPA ranging from Kalman filters [15]–[18], neural networks
[19], fuzzy logic [20], genetic algorithms [21], sliding mode observer [22], component adaptation
approach [23] - [24], and hybrid diagnosis [4]. Some of the linear approaches have been extended
to fully nonlinear models of the engine that have led to nonlinear extensions of the Kalman filters
to extended Kalman filters (EKF) and unscented Kalman filters (UKF) [12] and [25].
In this paper, a modular and a hierarchical non-interacting MM-based approach is proposed
for both single and concurrent permanent sensors fault detection and isolation in a single spool
jet engine. In our proposed hierarchical approach, multiple levels of detection filters are used to
detect and isolate concurrent faults where at each time only one level of the scheme is active
according to the engine sensors health status.
We also propose a hybrid Kalman filter (HKF) [17] scheme that consists of a single nonlinear
on-board engine model (OBEM) augmented with piecewise linear (PWL) models constituting
as the MM-based estimators to cover the entire engine operating regime. Also, the discrete-time
derivation of the HKF is formally shown in this paper. Therefore, multiple HKFs are constructed
to ensure that the FDII algorithm works effectively in a wide range of operating conditions by
decomposing the engine operating range into sub-regions each represented by a PWL model.
We then apply a Bayesian approach to generate a general combined model based on the PWL
models normalized weights. This provides us with a soft transition or interpolation among the
PWL models. At any given operating point, the combined model represents the nonlinear system
dynamics with less error than a single PWL model.
For interpolating the PWL models, various techniques are available in the literature including
Takagi-Sugeno fuzzy models [26], the functional state approach using discrete representations
such as automaton for describing the transitions between regions [27], Markov mixtures of
experts to define the transitions [28], etc. Unlike the functional state approach, the Bayesian
and fuzzy approaches enable soft transitions among the PWL models as the operating condition
changes between various points. Moreover, the Bayesian approach uses the statistical method to
infer the applicable and most appropriate operating regime at each time instant.
Sensor fault severity identification is another important part of our proposed FDII framework
which may either be utilized for reconstructing the correct sensor measurements to design an
output feedback control strategy or to remove the faulty sensor from the measurement set. The
reconstructed fault-free system outputs are essential for a normal closed-loop system operation
with a faulty sensor in order to maintain the system stability and set-point tracking under faulty
situations [29]. Towards this end, we have integrated the generalized likelihood ratio (GLR)based method [30]–[32] with our MM-based scheme to estimate the sensor fault severity under
various single and concurrent fault scenarios. The MM-based scheme does provide us with the
time and location of a sensor fault that subsequently simplifies the use of the GLR method
without requiring to solve an optimization problem for determining the fault detection time as
part of its fault estimation process.
August 28, 2015

DRAFT

3

One of the major non-fault related factors leading to false alarms in a given FDI scheme arises
due to the normal aging effects of the engine components such as compressor fouling, turbine
erosion and wear or buckling of a component. The above phenomena are associated with gradual
degradation of the health parameters from their healthy reference baselines, which is a type
of modeling uncertainty. Different approaches, including the robust generation and evaluation
methods, are introduced in [33] to improve the robustness of the on-line FDI algorithm with
respect to modeling uncertainties. In [34], modeling uncertainty is assumed as the structured input
with known distribution matrix and certain rank conditions in order to design a decoupling filter
to achieve the robustness. However, in most applications including ours, usually it is challenging
to model the uncertainty in a structured format. Therefore, robust residual evaluation methods in
which the adaptive threshold is designed based on the upper bound of the unstructured modeling
uncertainty that use the historical data have been introduced. Nevertheless, it is required to
design an individual threshold for each fault type. In addition, in few works such as [4], a
periodic updating mechanism is used for the on-board engine model in order to enhance the
robustness and the fault sensitivity.
In our paper, it is assumed that the reference baselines for the engine health parameters, namely
the efficiencies and the mass flow rates of the turbine and the compressor, can be estimated after
a certain number of flights by using an off-line health monitoring module [35] for determining
the percentages of degradations. The applied method can be utilized off-line in a ground station
given that the health parameters are degraded more slowly as compared to an abrupt fault.
Consequently, in order to ensure the reliability of our proposed FDII scheme during the engine
life cycle, the reference baselines are periodically updated for the OBEM health parameters.
This procedure will be useful in preventing occurrence of false alarms due to the engine health
parameters degradations. Despite this updating process, there is always a mismatch between the
values of the real engine health parameters and the ones that are used in the OBEM due to the
off-line estimation errors of the health monitoring method. Therefore, we will also investigate
the robustness of our proposed FDII scheme with respect to various estimation errors through
performing Monte Carlo simulations.
To summarize the main contributions of this paper can be stated as follows:
1) The explicit derivations of the discrete-time HKF scheme is formally provided in Appendix
A to complement the continuous-time HKF scheme derivations provided in [17].
2) A modular and hierarchical non-interacting MM-based HKF structure is developed to detect
and isolate single and concurrent sensor faults during the entire engine operating regime
(flight profile) having lower fault detection time and better robustness towards the engine
health parameters degradations as compared with the other linear and nonlinear filtering
methods in the literature. Our method is capable of operating during the entire engine life
cycle through periodically updating the reference baselines of the engine health parameters
for the nonlinear OBEM. In addition, the computational time of our proposed method is
lower than the investigated linear and nonlinear methods.
3) Through performing extensive simulation studies and measuring the mean of residual
signals, it is shown that the HKF is capable of estimating the engine outputs more
accurately with less number of operating points and false alarm rates as compared to
multiple linear Kalman filters (MLKF) based FDI method.
4) Our proposed MM-based HKF structure is integrated with the GLR scheme to estimate
the fault severity. This is accomplished by specifically eliminating the estimation of the
fault detection time which is a necessary step in the standard GLR scheme.
August 28, 2015

DRAFT

4

5) Our proposed multiple HKF-based FDI scheme is compared with the MM-based schemes
that utilize various linear and nonlinear filtering approaches such as the linear Kalman filter
(LKF), extended Kalman filter (EKF), unscented Kalman filter (UKF) and the cubature
Kalman filter (CKF) in terms of the promptness of the fault detection, false alarm and
incorrect fault detection rates, robustness with respect to the engine degradations and
computational time.
This paper is the extended version of [36] from both theoretical and case-based simulation
scenarios perspectives. Theoretically, we formally show how the HKF is derived and how
the previously proposed sensor FDI algorithm in [36] is integrated with the GLR method to
investigate the sensor fault identification problem. The case study scenarios are also extended
from the following perspectives, namely (a) different single and concurrent sensor faults are
considered, (b) the reductions in the fault detection times, lower false alarms and incorrect fault
detection rates under various parameter uncertainty levels are investigated, and (c) our proposed
FDI scheme is compared with linear and nonlinear observer-based methods in the literature
(LKF, EKF, UKF, and CKF) in terms of fault detection time, robustness towards the engine
health degradations and the computational time.
The remainder of the paper is organized as follows. Section II provides the design of the
HKF which includes the OBEM and the PWL models. The interpolation of the PWL models is
provided in Section III, and Section IV describes the overall hierarchical MM-based FDI scheme
applicable to both single and concurrent sensor faults. The modified GLR-based approach is
developed and described in Section V for estimating the severity of a sensor bias fault. A number
of Monte Carlo simulation scenarios and cases associated with various faults under multiple
operational conditions are presented. Finally, comparisons with various linear and nonlinear
filtering methods (LKF, EKF, UKF and CKF) are included in Section VI, and the paper is
concluded in Section VII.
II. HYBRID KALMAN FILTER (HKF) DESIGN
An actual aircraft gas turbine engine, used for deriving the on-board engine model (OBEM),
and which is used for on-line diagnostic analysis can be described according to the following
representation:
Ẋ(t) = Fc (X(t), H(t), U (t), ζ(t)),
Y (t) = Gc (X(t), H(t)) + v(t),

(1)

where X(t) ∈ Rn , H(t) ∈ Rr , Y (t) ∈ Rq , U (t) ∈ Rp , ζ(t) ∈ R2 and v(t) ∈ Rq denote
the engine state variables, health parameters, sensor measurements, input signals, and Gaussian
zero-mean process and measurement noise at time t, respectively. The engine dynamics is also
an implicit function of the ambient condition parameters including the ambient temperature and
pressure. The ambient parameters are defined in terms of the environmental parameters including
the altitude and the Mach number. The engine health parameters in (1) will become degraded
from their healthy reference baselines during the entire engine life cycle. Moreover, it is assumed
that the inputs and outputs of an actual engine are discretized with sufficiently small sampling
period for performing simulation and implementation of our proposed on-line fault detection,
isolation and identification (FDII) scheme. The discrete-time representation of the variables in
(1) are denoted by X(k), Y (k), H(k), U (k), ζ(k) and v(k) with the discrete dynamic functions
F and G replacing Fc and Gc , respectively.
August 28, 2015

DRAFT

5

The hybrid Kalman filter (HKF) consists of two main blocks that include a nonlinear on-board
engine model (OBEM) and multiple piecewise linear (PWL) models derived at different operating
points to cover the entire engine operating range. The continuous-time fault-free representation
of the OBEM which can also be derived based on thermodynamics laws is now given as follows:
ẊOBEM (t) = fc (XOBEM (t), HOBEM (t), U (t)),
YOBEM (t) = gc (XOBEM (t), HOBEM (t)),

(2)

where XOBEM (t) ∈ Rn and YOBEM (t) ∈ Rq denote the OBEM state variables and outputs. Both
the actual engine and the OBEM operate in parallel under the same flight conditions. Moreover,
HOBEM (t) ∈ Rr denotes the OBEM health parameters that can be represented by HOBEM (t) =
λT h(XOBEM (t)), in which λ defines the OBEM health parameters degradation factors (reference
baselines). This parameter is an all-ones vector for the OBEM that represents a healthy or nondegraded engine. It is also assumed that the OBEM health parameters can be periodically updated
and λ is considered as fixed in between the updating intervals. In addition, h(XOBEM (t)) denotes a
state-dependent smooth function that corresponds to the compressor and the turbine performance
maps and is modeled as a polynomial function for our developed gas turbine engine model. The
updating process for the OBEM health parameters will be described below in detail. Moreover,
for performing simulation and implementations, the continuous-time OBEM model is assumed
to be discretized with sufficiently small sampling period. The discrete-time representation of the
variables in (2) are denoted by XOBEM (k), YOBEM (k), HOBEM (k) and U (k) with the discrete
dynamic functions f and g replacing fc and gc , respectively.
The continuous-time OBEM model is linearized and discretized at multiple operating points
(corresponding to engine steady-state values) that are denoted by (Xssi , Ussi , Yssi ) with sufficiently
small sampling period. The constructed multiple linear discrete-time state-space models is now
given by:
∆Xi (k + 1) = A|Xssi ∆Xi (k) + B|Ussi ∆Ui (k),
∆Yi (k) = C|Xssi ∆Xi (k),

(3)

where i ∈ 1, . . . , L (L is the number of the operating points), A|Xssi , B|Ussi and C|Xssi denote the
state-space matrices associated with the ith operating point, and ∆Xi (k) = XOBEM (k) − Xssi ,
∆Yi (k) = YOBEM (k) − Yssi and ∆Ui (k) = U (k) − Ussi . In the Appendix A, the linearization and
discretization process of the continuous-time OBEM model is described in more detail. In this
paper, it is assumed that the OBEM linearized models are obtained with the health parameters that
are set to their healthy reference baselines. Also, the health parameter effects due to changes in
the engine state variables have been implicitly incorporated into the matrices A and C, although
the deviations from the healthy reference baselines are not incorporated in the linear models.
For each linearized model, an off-line linear Kalman filter is designed to estimate both the
actual engine states and sensor outputs as follows:
∆X̂i (k + 1) = A|Xssi ∆X̂i (k) + B|Ussi ∆Ui (k) + Kssi (Y (k) − Ŷi (k)),
∆Ŷi (k) = C|Xssi ∆X̂i (k),

(4)

where ∆X̂i (k) = X̂(k) − Xssi , ∆Ŷi (k) = Ŷ (k) − Yssi and Kssi denotes the steady-state Kalman
filter gain matrix. For the purpose of constructing the multiple model HKFs, the steady-state
Kalman gain matrices as well as the matrices A|Xssi and C|Xssi that are constructed associated
August 28, 2015

DRAFT

6

with multiple operating points are stored in a look-up table. According to (4), the linear Kalman
filter does not take into account the effects of the health parameter degradations from their healthy
reference baselines since the state-space matrices have already been determined for an all-ones
λ. Consequently, it does not have the required level of robustness for handling a vast number
of health degradations that occur during the entire engine life cycle. Therefore, it is essential
to update the health parameter reference baselines of the OBEM to maintain the reliability and
accuracy of the state estimates and the performance of the FDII scheme through out the entire
engine life cycle operation. For this purpose, the health parameter reference baselines can be
estimated by an off-line health monitoring system and then periodically updated in the discretetime OBEM model. Therefore, the on-line FDII scheme is integrated with an off-line health
monitoring system.
The frequency of the health monitoring system updates is significantly lower than that of
the on-line FDII algorithm, since the health degradation process is slow and gradual during
one flight, although their accumulated effects after a number of flights may generate a large
discrepancy between the OBEM and the actual engine. The off-line health monitoring module
can either be a single augmented Kalman filter that estimates the health parameters based on
the collected data during several flights or a nonlinear approximation method such as a neural
network that receives the engine historical data. It must be noted that both the health monitoring
and updating mechanism can be performed on-line having a sufficiently large sampling interval.
However, the health parameter estimation process needs to be terminated whenever a fault occurs
in the engine to avoid generating incorrect estimates of the health parameters.
The off-line health monitoring module uses the collected input and the measured output data
to estimate the health parameter reference baselines, that is λ̂. The estimated λ̂ is updated
periodically in the OBEM model within a determined time interval T that is usually a given
number of flights or days. The update process feeds λ̂ into the OBEM model so that the online FDII scheme can operate within the neighborhood of the degraded actual engine condition.
Consequently, the OBEM model can now be re-written as follows:
XOBEM (k + 1) = f (XOBEM (k), λ̂h(XOBEM (k)), U (k)),
YOBEM (k) = g(XOBEM (k), λ̂h(XOBEM (k))).

(5)

It must be noted that the OBEM model that is provided in (2) is used only once for performing
the linearization process with λ set to an all-ones vector, although the OBEM model in (5), with
its reference baselines periodically updated, is the one that is utilized in our proposed HKF
structure as well as in our on-line FDII scheme. For sake of notational simplicity, we use the
same notations for the above two versions of the OBEM model.
To construct the HKF representation we modify (4) where the steady-state variables are
replaced by the OBEM states and outputs that are obtained from (5) and also by using the
previously stored steady-state Kalman filter gain and state-space matrices as follows:
X̂i (k + 1) − XOBEM (k + 1) = A|Xssi (X̂i (k) − XOBEM (k)) + Kssi (Y (k) − Ŷi (k)),
Ŷi (k) = C|Xssi (X̂i (k) − XOBEM (k)) + YOBEM (k).

(6)

In the above model the effects of the input and the B matrix are eliminated from the HKF
formulation given the fact that these have already been accounted for by the OBEM model. The
procedure for derivation of the discrete-time HKF is formally shown in the Appendix A.
August 28, 2015

DRAFT

7

Consequently, multiple HKFs are designed by using (6) for multiple operating points. After
updating the OBEM health parameter reference baselines, there is no longer a need to recalculate
the A|Xssi , C|Xssi and Kssi matrices for each operating point, since the effects of the health
parameter degradations have been incorporated into the XOBEM and YOBEM as given by (5).
One of the factors that can affect the efficiency of our proposed FDII algorithm is the
estimation error of the off-line health monitoring module, that leads to mismatches between the
actual engine and the OBEM outputs. Larger estimation errors may increase the fault detection
time, and lead to occurrence of false alarms and incorrect fault detection rates. The acceptable
ranges of the estimation errors that do not lead to false alarms and incorrect fault detection
rates will be specified subsequently in Section VI for different health parameters under various
healthy and faulty scenarios. Moreover, the reliability of our proposed on-line FDII strategy as a
function of different mismatching factors between the actual operational engine and the OBEM
(as represented by the reference baseline estimation errors (RBEE = |λ−λ λ̂| )) and the magnitude
of the process and measurement noise signals will be investigated subsequently in Section VI
by means of a confusion matrix analysis.
In order to develop our proposed HKF-based scheme as an FDII strategy, multiple piecewise
linear (PWL) models need to be generated for each fault hypothesis at various operating points.
Moreover, the PWL models will be integrated and fused to cover the entire operational regime
of an engine. The detail description of this process is provided below.
III. P IECEWISE L INEAR M ODELS (PWL) I NTERPOLATION
Any given linear model of a nonlinear system has a limited operating range in which it remains
valid. Nevertheless, our ultimate goal is to obtain a globally valid model which is valid for the
entire operating regime. Therefore, the full operating range is divided into several sub-regions
where each is defined around an operating point for which a piecewise linear (PWL) model
can be derived [37]. The PWL models can then be integrated in order to construct a parametervarying general model whose parameters are the PWL models weights that are obtained through
an on-line Bayesian approach. This will provide one with a soft interpolation among the PWL
models as opposed to a hard switching among them. In this paper, the engine inputs including
the fuel flow rate and ambient variables are used to partition the engine operational regime into
multiple operating points for which the PWL models are constructed. These operating points are
associated with different flight conditions such as take-off, climbing, cruise and landing modes.
The selected number of operating points depends on (a) the required HKF state estimation
accuracy, (b) the FDII strategy reliability on correct decisions, and (c) false alarm rates within
the range of the applied health parameter degradations. Hence, if there are no concerns on the
memory utilization, the number of the operating points can be selected to be as high as possible
to enhance the HKF estimation accuracy and the FDII scheme valid decision rates and also to
decrease the false alarm rates.
One of the important advantages of our proposed HKF scheme is in requiring a smaller
number of operating points as compared to standard linear Kalman filters for covering an entire
operational regime of the engine. This is facilitated and made possible due to substitutions of
the steady-state variables in (6) by the OBEM state and output variables. The operating range
of a PWL model in (4) is only limited to the neighborhood of a corresponding operating point,
although this can be extended to a larger range in (6) given that XOBEM and YOBEM are changed
according to the engine operating condition, which enable the PWL model to be valid in a wider
range.
August 28, 2015

DRAFT

8

Since, the sensor fault is injected into the actual gas turbine engine and not the OBEM, we
have for the faulty engine
X(k + 1) = F(X(k), H(k), U (k), ζ(k)),
Y (k) = G(X(k), H(k)) +

q
X

bs zs δ(k − kf s ) + v(k),

(7)

s=1

where F and G represent the discrete-time dynamic equations of the actual gas turbine engine, q
is the number of sensors, bs represents the sth sensor bias fault magnitude and zs represents the
fault location vector that has a unit value for the sth element while the other elements are set to
zero, and δ(k − kf s ) denotes a unit step function that occurs at the sample kf s corresponding to
the sth fault occurrence time. The bias is set to zero for the healthy sensor scenario. Therefore,
there are a total of q + 1 sensor modes (corresponding to one healthy and q faulty sensor modes).
The PWL models constructed for multiple operating points are now used to compute the
corresponding Kssi for various sensor modes of the (7). The matrices A|Xssi , C|Xssi and Kssi are
finally stored in a look-up table and are used to construct multiple HKFs (MHKFs) as given by
(6) for all the sensor modes. Therefore, the HKF for the j th sensor mode at the ith operating
point is designed as follows:
X̂ (i,j) (k + 1) − XOBEM (k + 1) = Ai (X̂ (i,j) (k) − XOBEM (k)) + Kssi (Y (k) − Ŷ (i,j) (k)),
Ŷ (i,j) (k) = C i (X̂ (i,j) (k) − XOBEM (k)) + YOBEM (k) + bdj aj δj (k),

(8)

where i = 1, . . . , L, j = 1, . . . , (q + 1), bdj denotes the pre-determined sensor bias fault that can
be different from the actual sensor fault bs that is injected into (7), aj denotes the q-dimensional
vector and is one of q + 1 modes of a which is the fault parameter vector. For the healthy mode
or j = 1, the fault parameter vector, a, is set to a zero vector whereas for j = 2, . . . , q + 1, aj
has a unit value for the (j − 1)th element and all the other elements are set to zero. Moreover,
Ai = A|Xssi , C i = C|Xssi and Kssi are the previously stored state-space and Kalman gain matrices
that depend on the ith operating point. Therefore, L × (q + 1) HKFs are constructed covering the
entire engine operating range corresponding to different sensor modes. It must be noted that the
state-space matrices depend only on the operating points but the Kalman gain matrices depend
on both the operating points and the noise covariance matrices that are the same for all the
sensor fault modes, given that the OBEM does not take into account the effects of sensor faults.
Although the operating range of a PWL model in the HKF scheme has been increased as
compared to that of the one that uses the standard linear Kalman filters by replacing the steadystate values with the OBEM variables, still none of the PWL models are solely valid over the
entire operating range of an engine. Therefore, corresponding to each PWL model one can
associate a validity function that is based on its normalized weight as obtained by means of the
Bayes formula. For this purpose, the residual vectors γ (i,j) and the covariance matrices S (i,j)
that are generated by the multiple HKFs are used to compute the likelihood function f (i,j) for
the j th sensor mode at the ith operating region as follows:
γ (i,j) (k) = Y (k) − Ŷ (i,j) (k), S (i,j) (k) = cov(γ (i,j) (k)),
−1
1
p
× exp[ (γ (i,j) (k))T (S (i,j) (k))(−1) (γ (i,j) (k))],
f (i,j) (γ (i,j) (k)) =
2
(2π)q/2 |S (i,j) (k)|

August 28, 2015

(9)

DRAFT

9

where it is assumed that the innovation sequence generated by the hybrid Kalman filter, γ (i,j) (k),
is a Gaussian white noise process with zero mean and covariance matrix S (i,j) (k) that is calculated
numerically. The normalized weights for the j th sensor mode are updated recursively by using
the Bayes formula as follows:
f (i,j) (γ (i,j) (k))w(i,j) (k − 1)
.
w(i,j) (k) = PL
(i,j) (γ (i,j) (k))w (i,j) (k − 1)
i=1 f

(10)

The weights computed above should also remain outside a narrow bound to avoid becoming
close to zero as:
if w(i,j) (k) > ρ then w(i,j) (k) = w(i,j) (k),
if w(i,j) (k) ≤ ρ then w(i,j) (k) = ρ,

(11)

where ρ is a design parameter that is determined by trial and error and it invokes the PWL
models that have very small weights and it avoids them from being removed from the set of L
models. It is also useful for numerical robustness of the recursive weight algorithm [38].
Following the computation of the normalized weights, multiple time-varying models are now
constructed that cover the entire engine operating range subject to various sensor modes. For each
sensor mode, the corresponding model state-space matrices as well as the weighted innovation
vector and the covariance matrix are obtained by using the PWL model state-space matrices and
their associated normalized weights as follows:
Ajc (k)

=

L
X

w

(i,j)

i

(k)A ,

Ccj (k)

=

i=1

γcj (k) =

L
X
i=1

Ajc (k)

L
X

w(i,j) (k)C i ,

i=1

w(i,j) (k)γ (i,j) (k),

Scj (k) =

L
X

(w(i,j) (k))2 S (i,j) (k),

(12)

i=1

Ccj (k)
th

denote the weighted state-space matrices of a linear time-varying model
and
where
associated with the j sensor mode. Also, γcj (k) and Scj (k) denote the weighted innovation
vector and the covariance matrix of the j th sensor mode, respectively. The above procedure is
now designated as the PWL models interpolation. Consequently, (q + 1) weighted innovation
vectors and covariance matrices that operate through out the entire engine operational regime
are used in the next two sections to develop our proposed FDII scheme.
IV. SENSOR FDI VIA MULTIPLE-MODEL-BASED SCHEME
In this section, the overall structure of our proposed MM-based FDI scheme is presented.
It is assumed that the fault parameter vector can take on only one of (q + 1) sensor modes
as aj . Therefore, at each operating point there are (q + 1) PWL models; one for the healthy
sensors scenario and q corresponding to various faulty sensor scenarios, that have been designed
and integrated with only one OBEM for constructing the multiple HKF (MHKF)-based scheme
as formulated in (8). The innovation vectors and the covariance matrices that are generated
by the MHKFs are fused and weighted as given by (12). Finally, there are (q + 1) weighted
innovation vectors and covariance matrices that are used in the MM-based FDI scheme, where
they are operating under different healthy and faulty sensor scenarios during the entire engine
operational regime.
August 28, 2015

DRAFT

10

Y(k)
U(k)

OBEM
𝑿
++

𝒀𝑶𝑩𝑬𝑴 (𝒌)
(𝒌 + 𝟏)

(𝒊,𝒋)

𝑿
𝒁−𝟏

𝜸 𝟏,𝒋 (𝒌), 𝑺 𝟏,𝒋 (𝒌)
(𝒌)
+-

𝜸𝒄 𝟏 (𝒌), 𝑺𝒄 𝟏 (𝒌)

𝒘(𝟏,𝒋)
𝑳

++

𝑪𝒊

-+

𝜸 𝟐,𝒋 (𝒌), 𝑺 𝟐,𝒋 (𝒌)

𝒘 𝒊,𝒋

𝒘(𝟐,𝒋)

𝒊=𝟏

+

Ai

𝑲𝒔𝒔 𝒊
Offline PWL Kalman Filter

𝜸

𝒊,𝒋

(𝒌)

𝑺

𝒊,𝒋

PWL(i , j) Model from
Look - up Table with
dimension (L , q+1)

𝒊,𝒋
𝒊,𝒋
(𝒌) 𝜸 (𝒌), 𝑺 (𝒌)

𝒘(𝒊,𝒋)

𝜸 𝑳,𝒋 (𝒌), 𝑺 𝑳,𝒋 (𝒌)

𝒘(𝑳,𝒋)

+
++

Combination of
𝒌 𝜸 𝒊,𝒋 (𝒌)
L linear Models
for jth Sensor
Mode
𝒋

𝑳

(𝒘 𝒊,𝒋 (𝒌))𝟐 𝑺 𝒊,𝒋 (𝒌)
𝒊=𝟏

𝜸𝒄 𝟐 (𝒌), 𝑺𝒄 𝟐 (𝒌)
𝜸𝒄 𝒋 (𝒌), 𝑺𝒄 𝒋 (𝒌)

𝒋

𝜸𝒄 𝒌 , 𝑺𝒄 (𝒌)

𝑷𝟏
𝑷

𝟐
Conditional
Hypothesis
Probability 𝑷
𝒋
Calculator
Equations
(13) and (14) 𝑷𝒒+𝟏

MAX

𝝀

With Updating
Interval T

(𝒊,𝒋)

Bayesian Weight Estimation

𝑿𝑶𝑩𝑬𝑴 (𝒌)

Fault
Detection
&
Isolation

𝜸𝒄 (𝒒+𝟏) (𝒌), 𝑺𝒄 (𝒒+𝟏) (𝒌)

Fig. 1: The MM-based FDI scheme based on the HKF corresponding to multiple operating
regime.
A. Single Fault Detection and Isolation (FDI) Scheme
In the MM-based approach [10]–[12], the hypothesis conditional probability Pj (k) is defined
as the probability that the fault parameter a assumes the mode aj , j = 1, . . . , q+1, conditioned on
the observed measurement history up to the k th sample, that is: Pj (k) = Pr[a = aj |Y(k) = Yk ],
where Y(k) is the random vector measurement history with Y (1), Y (2), . . . , Y (k) partitions displaying the available measurements up to the k th sample time. Similarly, Yk is the measurement
history vector realization that has the partitions of Y1 , Y2 , . . . , Yk . Therefore, the conditional
probability can be computed recursively as follows:
fY (k)|a,Y(k−1) (Yj |aj , Yk−1 )Pj (k − 1)
,
Pj (k) = Pq+1
~=1 fY (k)|a,Y(k−1) (Yj |a~ , Yk−1 )P~ (k − 1)

(13)

where fY (k)|a,Y(k−1) (Yj |aj , Yk−1 ) denotes the Gaussian density function for the current measurement given by:
fY (k)|a,Y(k−1) (Yj |aj , Yk−1 ) =
(2π)q/2

1
q

× exp[
Scj (k)

−1 j
(γ (k))T (Scj (k))(−1) (γcj (k))],
2 c

(14)

where γcj (k) and Scj (k) are given by (12).
If the j th sensor mode occurs, the probability associated with the j th model will be larger
than that of the others since its corresponding innovation vector and also the determinant of the
covariance matrix will be much smaller than those that are predicted by the other filters and
which are mismatched with the assumed fault scenario. Hence, the condition of the system and
the location of a single faulty sensor can be detected and isolated based on evaluating Pj (k) and
determining its maximum value. Consequently, our proposed MM-based approach is capable of
detecting and isolating sensor faults. Fig. 1 shows the structure of the MM-based FDI scheme
that employs the MHKFs for the entire engine operational regime.
Our contribution here is the modification of our previously developed MM-based structure
in [10] and [11] that utilizes only multiple standard linear Kalman filters and is designed for a
single operating point into a general strategy that is applicable to the entire engine operational
regime by means of the multiple HKF scheme and fusion and integration of the corresponding
innovation vectors and covariance matrices.
As we will describe in more detail in Section VI, in this paper five sensor faults are considered.
Therefore, the total number of modes corresponding to each operating point of the jet engine
is six, where mode #1 (P1 ) corresponds to the healthy engine sensors and modes #2 to #6 (P2
August 28, 2015

DRAFT

11

TABLE I: Operational modes corresponding to various possible two concurrent sensor faults.
Level
First

#1
Healthy
P2
P3

Second

P4
P5
P6

Operational Modes
#2 #3 #4
P2
P3
P4
P2
P2
P2
P2
P3
P4
P3
P3
P3
P2
P3
P4
P4
P4
P4
P2
P3
P4
P5
P5
P5
P2
P3
P4
P6
P6
P6
P2
P3
P4

#5
P5
P2
P5
P3
P5
P4
P5
P5
P5
P6
P5

#6
P6
P2
P6
P3
P6
P4
P6
P5
P6
P6
P6

to P6 ) correspond to a 3% sensor bias fault injected into sensors measuring the compressor exit
temperature (P2 ) and pressure (P3 ), shaft rotational speed (P4 ), the turbine exit temperature (P5 )
and pressure (P6 ), respectively.
B. Concurrent Fault Detection and Isolation (FDI) Scheme
For detection and isolation of two concurrent faults in the gas turbine engine, a hierarchical
MM approach is proposed in [9] and [11] as shown in TABLE I. In this scheme, it is assumed that
the engine starts operating with healthy sensors when the first level of filters are active and the
FDI scheme observes the sensors condition for occurrence of one of the five faulty modes. The
active bank of filters operate based on the weighted innovation vectors and covariance matrices
that are obtained in Section III.
Once the first fault is detected and isolated according to the maximum probability criterion,
the FDI scheme will activate the second level of the filters as shown in TABLE I for detection
of the second concurrent faults. In our hierarchical MM scheme, it is assumed that the sensor
faults do not occur simultaneously and there exists always a minimum time interval between
the occurrence of two sensor faults that are called concurrent faults. TABLE I indicates all
the possible configurations for the second bank of filters that are considered in the simulation
results in Section VI. The first filter in the second level always corresponds to the detected faulty
mode in the first level. For example, if the first filter detects a 3% bias in the compressor exit
temperature sensor (P2 ), then the first and the second filters in the second level respectively
correspond to P2 and the possibility of a larger bias fault within the compressor exit temperature
sensor. Moreover, the third filter corresponds to the concurrent bias faults in the compressor exit
temperature and pressure sensors (P2 and P3 ), the forth filter corresponds to P2 and P4 , etc. This
structure can easily be extended to the third and higher levels that correspond to occurrence of
multiple concurrent sensor faults. It must be noted that when the new bank of filters is activated
in the second level, the first bank will be disabled in order to avoid adding any unnecessary
computational burden. Therefore, at any given time only six filters are operating on-line.
Remark 1. In this paper, the value of the pre-determined sensor bias fault severity, bdj as
given by (8), is considered to be the same for all the filters and only one level of fault severity is
considered for the hierarchical MM scheme, although various pre-determined sensor bias faults
can easily be incorporated into our strategy by correspondingly increasing the number of filters
that are used in this scheme.

August 28, 2015

DRAFT

12

V. SENSOR FAULT IDENTIFICATION/ESTIMATION
One of the important requirements of a general control system is the capability to integrate the
FDI scheme with a fault identification or estimation module in order to estimate the severity of
a fault that has occurred in different components of the system such as sensors and actuators. In
this paper, a modified version of the generalized likelihood ratio (GLR) scheme is developed to
estimate the severity of a sensor bias fault. The GLR was initially proposed by Willsky in [30]
and was subsequently modified in [31], [32], [39] and [40]. This method is capable of estimating
the time, location and severity of an occurred fault using a selected threshold. However, in this
paper the occurrence time and location of a fault have already been determined by means of our
proposed MHKF-based approach as presented in Section IV. Hence, the fault severity will be
estimated through the development of a modified version of the GLR scheme. The GLR method
is a detection-estimation scheme in which the time of a fault is an a priori information input
to the estimation stage. In this paper, the GLR method is simplified by removing the detection
process through integrating the GLR method with our proposed MM-based FDI scheme, and
hence there is no longer a need to select a threshold for our FDII scheme.
Let us assume that a bias fault has occurred in the sth sensor with the severity bs at the sample
kf s and it is detected and isolated at the sample kds when the P1 (k) mode probability intersects
with the P(s+1) (k) mode probability. Also, for k ≥ kds , P(s+1) (k) has the highest value among
the other mode probabilities as shown mathematically by (s + 1) = argj max Pj (k). Therefore,
the effects of the sth sensor fault still remain in the residual of the MHKF that is associated
with the healthy mode. There are L PWL models corresponding to the healthy mode that will
be used subsequently to estimate the fault severity as follows:
i
X̂ (i,1) (k + 1) − XOBEM (k + 1) = Ai (X̂ (i,1) (k) − XOBEM (k)) + Kss
(Y (k) − Ŷ (i,1) (k)),

Ŷ (i,1) (k) = C i (X̂ (i,1) (k) − XOBEM (k)) + YOBEM (k),

(15)

in which the j th index and bdj in (8) are equal to one and zero, respectively, that are designated
to the healthy sensor mode. Therefore, the residual vector of the MHKF that is designed for the
healthy sensor mode at any subsequent time can be expressed as follows:
γ (i,1) (k) = bs Gis (k, kds )zs + v(k)

for k ≥ kds ,

(16)

where zs is defined in (7) and Gis (k, kds ) denotes the failure signature matrix that provides one
with the information on the failure propagation through the filter, and which also depends on
both kds and the sample time k at which the set of L innovation vectors are computed for the
MHKFs associated with the healthy mode. The signature matrices for a sensor bias fault can be
computed by using the recursive relations for the healthy mode MHKFs as follows [39]:
Gis (k, kds ) = I − C i Ai Jsi (k − 1, kds ),
i
Jsi (k, kds ) = Ai Jsi (k − 1, kds ) + Kss
Gis (k, kds ),

(17)

where Jsi (k, kds ) denotes the signature matrix for the state correction. The signature matrices
Jsi (k, kds ) and Gis (k, kds ) are defined separately for the L operating points over a window of
duration [kds , kds + N ], where N denotes the data samples window length. Also, it is assumed
that Jsi (kds − 1, kds ) = 0. The weighted signature matrices for a detected fault that cover the
entire operating range of the gas turbine engine can be computed as follows:
August 28, 2015

DRAFT

13

Gs (k, kds ) =

L
X

w(i,1) (k)Gis (k, kds ),

i=1

Js (k, kds ) =

L
X

w(i,1) (k)Jsi (k, kds ),

(18)

i=1

where w(i,1) (k) denotes the assigned weight to the ith PWL model that is designed for the
healthy sensor mode. Therefore, the modified GLR test which consists of a maximum likelihood
estimation (MLE) of bs when kds is known is used. For performing the maximization process,
normally the log likelihood function is chosen as the GLR criterion as follows:
J =

kds
+N
X

T
γc1 (k)(Sc1 (k))(−1) γc1 (k)

−

kds
+N
X

[γc1 (k) − bs Gs (k, kds )zs ]T (Sc1 (k))(−1) [γc1 (k)

k=kds

k=kds

− bs Gs (k, kds )zs ],

(19)

where γc1 (k) and Sc1 (k) denote the weighted innovation vector and covariance matrix for the
healthy sensor mode that are computed from (12). The maximization of the log likelihood
function is associated with the minimization of the second summation in (19). It can easily be
shown that the optimal solution to the MLE problem above is the estimated fault b̂s = dcss , that
is the unbiased estimate of bs having the minimum variance as follows:
ds =

zsT

kds
+N
X

1
(−1) 1
GT
γc (k),
s (k, kds )(Sc (k))

k=kds

cs =

zsT

kds
+N
X

1
(−1)
GT
Gs (k, kds )zs .
s (k, kds )(Sc (k))

(20)

k=kds

The above procedure can also be employed for estimating the severities of multiple concurrent
faults that have already been detected and isolated by means of our proposed hierarchical MMbased FDI approach.
Remark 2. In this paper, it is assumed that there is no feedback and information sent from
the fault identification module to the FDI scheme; otherwise, the FDI scheme needs to have
a variable structure for updating the pre-determined sensor bias faults based on the estimated
faults severities.
VI. SIMULATION RESULTS
In this section, simulation results and performance evaluation of our proposed sensor FDII
scheme as applied to several fault scenarios are presented for a nonlinear mathematical model of
a commercial single spool jet engine previously developed in [12]. This model is generated based
on the rotor, volume and heat transfer dynamic behavior and is also validated by the commercially
available standard software GSP 10 [41]. The responses corresponding to our mathematical model
and the GSP match each other within an acceptable error tolerance (below 5%). For this study,
our model is extended to the entire flight profile and is simulated in the SIMULINK for use
as both the actual engine and the OBEM model. The actual engine operates at a given health
condition subject to the effects of process and measurement noise signals, whereas the OBEM
health parameters reference baselines are periodically updated to their recently estimated values
August 28, 2015

DRAFT

14

that are assumed to be generated by an off-line health monitoring module as described in Section
II.
Both the actual engine and the OBEM models operate in almost the same ambient conditions.
The set of nonlinear state-space equations that are used for the engine model are given as follows:
PCC
[(cp TC λṁC ṁC + ηCC Hu ṁf − cp TCC λṁT ṁT ) − cv TCC (λṁC ṁC + ṁf
TCC cv ṁCC
γRTCC
− λṁT ṁT )] +
(λṁC ṁC + ṁf − λṁT ṁT ),
VCC
ηm λṁT ṁT cp (TCC − TT ) − λṁC ṁC cp (TC − Td )
,
Ṅ =
π 2
JN ( 30
)
1
[(cp TC λṁC ṁC + ηCC Hu ṁf − cp TCC λṁT ṁT ) − cv TCC (λṁC ṁC + ṁf − λṁT ṁT )],
ṪCC =
cv mCC
RTM
β
ṖT =
(λṁT ṁT +
λṁ ṁC − ṁn ),
(21)
VM
β+1 C
ṖCC =

where X = [PCC , N, TCC , PT ]T denotes the state variable vector that includes the combustion
chamber pressure PCC and temperature TCC , the rotational speed N and the turbine pressure PT . There is a single actuator that supplies the fuel flow
well asi five
h (u = ṁf ) asγ−1
PCC
1
T
sensors measuring Y = [TC , PC , N, TT , PT ] , where TC = Td 1 + λη ηC [( Pd ) γ − 1] and
C
h
i
γ−1
T
TT = TCC 1 − ληT ηT [1 − ( PPCC
) γ ] denote the compressor and the turbine temperatures,
respectively. The variables PC and PT denote the pressures of the compressor and turbine,
respectively. Moreover, H = [ηC , ηT , ṁC , ṁT ]T denotes the health parameter vector, where
ηC and ηT denote the compressor and the turbine efficiencies, and ṁC and ṁT denote their
mass flow rates, respectively. The health parameters are also multiplied by their corresponding
health parameter degradation factors or the reference baselines λ = [ληC , ληT , λṁC , λṁT ]T .
The sensors are affected by the Gaussian measurement noise with the standard deviations of
SDv = [0.23, 0.164, 0.051, 0.097, 0.164]T corresponding to the percentages of output vector at
the cruise condition as defined in [12]. More details on other variables and constants in (21) are
provided in [12].
The flight condition is defined by two environmental variables, namely the altitude and the
Mach number. The ambient temperature and pressure can be computed according to Tamb =
Alt
Ts − 6.5Alt
and Pamb = Ps exp( −gM
), where Ts = 288 ◦ K and Ps = 1.01325 bar are set to the
1000
288R
standard condition, R = 8.31447 is the gas constant, g is the gravitational acceleration, and M
and Alt denote the Mach number and the altitude, respectively. Moreover, the ambient variables
are affected by the Gaussian process noise with the standard deviations of SDζ = [0.01, 0.01]T
corresponding to the percentages of standard conditions of ambient variables. The noise free
values of the ambient parameters are applied to the OBEM model. Also, the same control input
is applied to both the actual engine and the OBEM models. The system is simulated for 520
sec with the sampling rate of 0.01 sec. The profiles of the altitude, Mach number and the fuel
flow rate are shown in Fig. 2.
According to the overall structure of our proposed scheme, as shown in Fig. 1, six PWL
models are constructed (q = 5) corresponding to each operating point for the five faulty sensor
modes as well as the one healthy mode with the pre-determined sensor bias faults set to 3%
of the engine steady-state outputs under the cruise condition. The states and sensor outputs
August 28, 2015

DRAFT

15

2

0.3
0.25
0.2
0.15

x 10

4

1

1.8

0.9

1.6

0.8

1.4

0.7

Mach Number

0.35

Altitude (ft)

Fuel Flow Rate (Kg/m2)

0.4

1.2
1
0.8
0.6

0.1

0.4

0.05
0
0

200

300

Time (Sec)

(a)

400

500

600

0
0

0.5
0.4
0.3
0.2

0.2

100

0.6

0.1
100

200

300

Time (Sec)

400

500

600

0
0

100

200

300

400

500

600

Time (Sec)

(c)

(b)

Fig. 2: Profiles of (a) the fuel flow rate, (b) altitude, and (c) the Mach number during a flight
mission.
TABLE II: The operating point specifications corresponding to the designed PWL models, where
ṁf and Alt are measured in Kg/m2 and ft, respectively.
Models Corresponding
to Flight Conditions
Operating Point 1 (Climbing)
Operating Point 2 (Climbing)
Operating Point 3 (Cruise)
Operating Point 4 (Landing)
Operating Point 5 (Landing)

ṁf

M

Alt

0.38
0.38
0.25
0.3
0.3

0.2109
0.6585
0.85
0.5402
0.1203

4070.538
12708.33
16404.2
10424.87
2322.835

are estimated by means of our proposed MHKFs. Finally, six weighted matrices are calculated
according to (12) and are applied in the MM-based scheme for detecting and isolating the sensor
bias faults that occur at different points of the flight profile.
For our simulations, the measurement and process noise covariance matrices are set to 0.01I
and 0.1I, respectively, where I denotes an identity matrix. In this paper, the FDII scheme
is implemented for the entire flight profile including the climbing, cruise and landing modes.
Therefore, we assume and consider five operating points handling the entire flight profile (L = 5).
This is the minimum number of the operating points that are obtained for each faulty mode that
prevents the occurrence of a false alarm in the range of the applied health degradations. The
number of the operating points depends on the dynamics of the engine as well as the range of
the applied fuel flow rate and the environmental parameters. We simulated our proposed FDII
scheme with different number of operating points and also derived a confusion matrix for each
case in order to analyze the false alarm rates and to decide on the minimum required number of
the operating points that does not lead to any false alarms for the injected health degradations.
These details are not included here due to space limitations.
TABLE II shows the corresponding fuel flow rates as well as the flight conditions corresponding to all the applied PWL models that are designed for each sensor mode. In order to
efficiently track the variations of the system input during the climbing and landing conditions,
two operating points are selected for each of these two modes; whereas only one operating point
is selected for the cruise condition given the presence of a constant input.
During the fault detection and isolation process, a mode probability Pj (k) is generated for
each weighted model using (13) and by determining the maximum Pj (k), j = 1, . . . , q + 1, the
sensor fault is detected and isolated.
August 28, 2015

DRAFT

16

TABLE III: Fault detection time (FDT) corresponding to the maximum tolerable percentage of
the RBEE (∆λ) at different stages of the flight profile.
Faulty Scenario
Fault on TC
Fault on PC
Fault on N
Fault on TT
Fault on PT

Maximum Tolerable RBEE%
Compressor
Turbine
Compressor
Turbine
Compressor
Turbine
Compressor
Turbine
Compressor
Turbine

∆λṁC
∆λṁT
∆λṁC
∆λṁT
∆λṁC
∆λṁT
∆λṁC
∆λṁT
∆λṁC
∆λṁT

= 2.52
= 2.52
= 1.01
= 1.01
= 1.21
= 1.01
= 0.2
= 0.1
= 0.3
= 0.3

∆ληC
∆ληT
∆ληC
∆ληT
∆ληC
∆ληT
∆ληC
∆ληT
∆ληC
∆ληT

= 2.52
= 1.51
= 1.01
= 0.30
= 2.02
= 0.81
= 0.2
= 0.1
= 0.3
= 0.1

kf = 50
3.7
3.5
7.8
3.4
6.5
6.5
7.8
6
7
8

Sensor FDT (sec)
kf = 250
4.1
5.1
5.9
2.7
2.6
2.5
3
2.2
2.6
2.2

kf = 450
6.3
5.9
8
2.5
3.5
2.9
4.7
2.3
3.4
2.2

As described in Section I, the health parameter degradations due to the engine aging is one
of the non-fault related factors that should be considered for a reliable and an accurate FDII
strategy. Therefore, the performance of the FDII strategy is evaluated in presence of compressor
or turbine health parameters degradations from their healthy reference baselines. Generally, there
is a difference between the actual and the estimated health parameter reference baselines due to
off-line health monitoring estimation errors and also due to effects of noise and disturbances.
Therefore, the robustness of our proposed FDII scheme will also be evaluated subsequently in
this section with respect to the percentage of the reference baselines estimation errors.
A. Case 1: False Alarms Evaluation
Many factors such as (i) the dynamic mismatch between the OBEM model and the actual
engine model, (ii) the large estimation errors for an off-line health monitoring module, and
(iii) the process and measurement noise, may lead to false alarm flags. In order to evaluate the
reliability and the efficiency of our proposed FDII scheme in terms of avoiding false alarms, our
proposed scheme is simulated over the entire flight profile as shown in Fig. 2. The estimated
reference baselines (λ̂) may differ form the actual degraded engine (λ) and this can play a
source of uncertainty in our FDII scheme in addition to the process and measurement noise. It
should be noted that the dynamic mismatch between the actual engine and the OBEM models
is also another source of uncertainty, but this is not investigated in this paper. Therefore, the
robustness of the FDII scheme against the percentage of the reference baselines estimation error
RBEE (RBEE% = 100 × ∆λ
, where ∆λ = |λ − λ̂| and λ ∈ [0, 1]) is investigated in this section.
λ
Based on our simulation studies for the healthy sensors scenario, the FDII algorithm is robust
to the maximum percentage of reference baselines estimation error of 3% for the compressor
health parameters and a maximum of 2% for the turbine health parameters. Our proposed FDII
algorithm has declared no false alarms in the range of the above RBEEs, however the possibility
of false alarms will increase if the RBEEs are increased beyond these upper limits. In case that
one is within these limits, the mode probability of the healthy sensors is near one whereas the
other probabilities corresponding to the faulty sensor scenarios are all almost zero during the
entire flight profile. The values for the probabilities depend on the design parameter ρ. In general
the MM-based FDI is independent of ρ as long as it is selected as sufficiently small.
August 28, 2015

DRAFT

17

P1

P2

P3

P4

P5

P6

1

Mode Probability

Mode Probability

1

0.5

0
0

100

200

300

400

Time (Sec)

(a)

500

600

P

1

P

2

P

3

P

4

P

5

1.2

P

P

6

0.5

0
0

1

P

2

P

P

3

4

P

5

P

6

1

Mode Probability

1.2

1.2

100

200
300
400
Time (Sec)

(b)

500

600

0.5

0
0

100

200
300
400
Time (Sec)

500

600

(c)

Fig. 3: Mode probabilities for 3% bias fault applied at (a) kf = 50 sec to the TC sensor, (b)
kf = 250 sec to the PC sensor, and (c) kf = 450 sec to the N sensor in presence of the obtained
maximum tolerable percentage of the estimation errors in TABLE III for the compressor health
parameter reference baselines.

B. Case 2: A 3% Sensor Bias Fault Detection and Isolation
In this section, the performance of the sensor FDI scheme with respect to the fault detection
time and the robustness of the algorithm with respect to the percentage of reference baselines
estimation error (RBEE) are evaluated. This evaluation is performed during the entire flight
profile which lasts for 520 sec, when a pre-determined bias fault with the severity of 3% of the
engine steady-state output values under cruise condition occurs for a single sensor. TABLE III
shows the fault detection times for each single sensor fault scenario at different stages of the flight
profile including the climbing, cruise and landing modes. It also shows the maximum tolerable
percentages of the RBEE of engine health parameters that is denoted by ∆λ. It is assumed that
λ = 0.99. The limits are obtained through preserving two performance requirements for our
proposed FDI scheme, namely: 1) the sensor fault detection time (FDT) should be less than 8
sec, and 2) no false alarms should be generated. This implies that our proposed FDI scheme can
tolerate the reference baselines estimation error as long as the above performance requirements
are fulfilled.
In addition, Fig. 3 depicts the mode probabilities for three selected fault scenarios when the
bias fault occurs at different instants of the flight profile. The value of the injected fault is set
to 3% of the engine steady-state output values, while the percentage of the estimation errors for
the compressor health parameters reference baselines (RBEE) are set at their maximum tolerable
levels that are indicated in TABLE III.
By comparing the results in TABLE III, it can be concluded that the sensor fault detection
times during the cruise mode are much less than that of the other flight modes since there is
less variation of thrust and ambient conditions. In spite of a large input and ambient condition
variations during the climbing and the landing modes, it is still possible to detect a sensor fault
by applying our proposed FDI scheme. In order to show this capability, the sensor faults occur at
kf = 50 sec during the climbing mode, at kf = 250 sec during the cruise mode and at kf = 450
sec during the landing mode. The fault detection times are all indicated in TABLE III for various
fault scenarios.
To quantify the effectiveness and reliability of our proposed FDI scheme in presence of the
mismatch between the OBEM model and the actual engine, confusion matrices are obtained for
both the healthy and the faulty scenarios subject to different health parameter degradations and

August 28, 2015

DRAFT

18

TABLE IV: The confusion matrix for 3% RBEE of the compressor health parameters.
TC
PC
N
TT
PT
No Fault

TC
50
1
0
12
17
0

PC
0
49
0
0
0
0

N
0
0
50
0
0
0

TT
0
0
0
27
0
0

PT
0
0
0
10
29
0

No Fault
0
0
0
1
4
50

TABLE V: The confusion matrix for 3% RBEE of the turbine health parameters.
TC
PC
N
TT
PT
No Fault

TC
50
0
1
2
0
0

PC
0
50
0
0
1
1

N
0
0
48
0
0
0

TT
0
0
0
41
0
0

PT
0
0
0
0
44
0

No Fault
0
0
1
7
5
49

measurement noise. The rows in the confusion matrix are the fault conditions and the columns
are the actual isolated faults. The element in the ith row and j th column (CMij ) shows the rate
that fault j is isolated when fault i occurs. Ideally, the confusion matrix should be a purely
diagonal matrix. To obtain the confusion matrix for our sensor FDI scheme, 50 Monte Carlo
simulations are performed in which the percentage of the RBEE exceeds the maximum tolerable
limit that is reported in TABLE III. Tables IV and V show the confusion matrices for the healthy
and 3% faulty scenarios as applied to the gas turbine engine sensors during the cruise mode
subjected to 3% RBEE of both the compressor and the turbine health parameters, respectively.
In addition to the percentage of the RBEE, the measurement noise is another factor that
can affect the performance of our proposed FDI scheme. Consequently, TABLE VI depicts the
confusion matrix that is obtained for investigating the effects of measurement noise signals,
when their original standard deviations (SDv ) are multiplied by a factor of 20.
According to our simulations, the initial impact of any increase in either the standard deviation
of the noise measurements or the percentage of the RBEE is to delay the fault detection time.
Moreover, the occurrence of false alarms and incorrect fault detections are the other consequences
due to increases in the above uncertainty sources. Different performance indices can be defined
to quantify the robustness of our proposed sensor FDI algorithm with respect to the levels of
uncertainty sources [42], namely:
P5
P5 P5
P6
j=1
i=1 CMij (i 6= j)
j=1 CM6j
i=1 CMii
FPR = P6
, ACC = P6 P6
, IFDR =
, (22)
P6 P5
CM
CM
CM
6j
ij
ij
j=1
j=1
i=1
j=1
i=1
where FPR, ACC and IFDR denote the false positive (false alarm) rate, accuracy and incorrect
fault detection rate, respectively, that are calculated in TABLE VII to investigate the effects of
various uncertainty sources.
According to TABLE VII, increasing the levels of the uncertainty sources results in decrease
of ACC, but increase in FPR and IFDR. Notwithstanding these observations, our proposed sensor
FDI scheme still works sufficiently robust against a high level of measurement noise as well
as discrepancies between the OBEM model and the actual engine health parameters reference
August 28, 2015

DRAFT

19

TABLE VI: The confusion matrix for the measurement noise signals where the SDs are multiplied
by a factor of 20.
TC
PC
N
TT
PT
No Fault

TC
48
0
0
0
0
0

PC
0
46
0
0
0
0

N
0
0
49
0
0
0

TT
0
0
0
48
1
0

PT
0
0
0
0
47
1

No Fault
2
4
1
2
2
49

TABLE VII: Sensor FDI algorithm performance indices corresponding to different levels of
uncertainties.
Scenarios
3% Compressor RBEE
4% Compressor RBEE
3% Turbine RBEE
4% Turbine RBEE
Noise SD×20
Noise SD×25

FPR
0
0.12
0.02
0.04
0.02
0.02

ACC
0.85
0.66
0.94
0.91
0.956
0.91

IFDR
0.16
0.376
0.016
0.076
0.0004
0.012

baselines.
C. Case 3: Sensor Fault Detection and Isolation for Different Fault Severities
In real applications, there is no guarantee that the sensor bias fault severity always matches
the 3% pre-determined fault for which the MM structure is designed. Therefore, it is essential
to investigate the performance of our proposed MM-based FDI scheme for an applied sensor
bias fault having different severities starting from the minimum detectable bias. TABLE VIII
shows the average detection times for all the faulty modes as a function of the fault severities,
when the applied faults occur at different stages of the flight profile with also the maximum
tolerable RBEE% for the compressor health parameters that are indicated in TABLE III. It can be
observed from TABLE VIII that the higher the fault severity with respect to the pre-determined
3% fault, the later the detection time, given that the actual fault becomes further different from
the pre-determined bias fault. Note that the minimum detectable sensor bias fault is 2% that
requires larger time to be detected as compared to the higher fault severities, especially during
the climbing and landing flight modes.
It should be pointed out that our proposed MM-based structure is designed for a 3% sensor
bias fault and is not capable of detecting and isolating sensor bias faults that have far greater

TABLE VIII: The average sensor fault detection times for all fault modes as a function of the
fault severity at different stages of the flight profile, under the maximum tolerable RBEE% for
the compressor health parameters.
Fault Time
kf = 50 sec
kf = 250 sec
kf = 450 sec
August 28, 2015

2%
14.18
8.34
17.66

3%
6.56
3.64
5.18

4%
6.74
3.98
5.32

5%
7.42
4.76
6.66

6%
8.6
5.7
7.52
DRAFT

20

TABLE IX: The average of the maximum detectable sensor bias fault (as percentage of the
engine steady-state outputs) over the entire flight profile by using our FDI scheme when the
corresponding maximum tolerable RBEE% reported in TABLE III are applied to the compressor
and turbine health parameters.
Fault Modes
TC
PC
N
TT
PT

Max. Tolerable
Compressor RBEE%
25%
20%
12%
10%
70%

Max. Tolerable
Turbine RBEE%
70%
50%
100%
14%
75%

severities. For this purpose, TABLE IX shows the maximum detectable sensor bias faults by
using our proposed FDI scheme by also considering the maximum tolerable RBEE% for the
compressor and the turbine health parameters. The reported detectable fault severities in TABLE
IX correspond to the average fault severities for all the flight modes. However, if an applied sensor
bias fault increases beyond the corresponding maximum detectable fault as indicated in TABLE
IX, an incorrect fault may be detected. Therefore, it is recommended that one incorporates more
models within the MM-based structure corresponding to higher pre-determined sensor faults to
become capable of detecting and isolating faults with higher magnitudes in shorter durations of
time.
D. Case 4: Sensor Fault Severity Identification
One of the advantages of our proposed MM-based FDI scheme is its capability in providing
information on the occurred fault severity. However, this method is not precise for fault severity
estimation, since it provides one with only information on the level of an occurred fault instead of
the exact fault severity. Therefore, the modified GLR method introduced in Section V is integrated
with our MM-based FDI scheme to estimate the sensor bias fault magnitude and severity. In
order to investigate the performance of our proposed fault severity estimation method, the notion
of the percentage of weighted mean square normalized residual errors (WMSNE) are obtained
for various detected fault scenarios as follows:
WMSNE% =
PK
(i,s+1)
(i,s+1)
L
(k)( (Y (k)−YŶ (k) (k)) )2
1 X k=kds w
Pκ
× 100,
(i,s+1) (k)
L i=1
k=kds w
Ŷ (i,s+1) (k) = C i (X̂ (i,s+1) (k) − XOBEM (k)) + YOBEM (k) + b̂s zs δ(k − kds ),

s = 1, . . . , q (23)

where K is the total simulation samples and Ŷ (i,s+1) (k) is the numerically constructed output
associated with the (s + 1)th operational mode corresponding to the sth fault, in which the fault
vector effect is generated based on the estimated fault magnitude b̂s and the fault occurrence
detected time kds . It must be noted that there is no feedback and information sent from the fault
estimation module to the filters for updating their pre-determined bias faults. TABLE X shows
the average percentage of the WMSNE that is measured for different stages of the flight profile
as a function of the bias fault severity as applied to various sensors.
August 28, 2015

DRAFT

21

TABLE X: The average WMSNE% for all flight modes as a function of the bias fault severity
applied to various sensors without considering any estimation error for the health parameter
reference baselines.
Fault
TC
PC
N
TT
PT

2%
0.0085
0.1360
0.0111
0.0137
0.1313

3%
0.0067
0.1106
0.0109
0.0066
0.1269

4%
0.0073
0.1293
0.0107
0.0128
0.1359

5%
0.0124
0.1689
0.0105
0.0312
0.1503

6%
0.0195
0.2537
0.0104
0.0944
0.1615

According to TABLE X, the average percentage of WMSNE is less than 0.5% for various
faults that are limited within the operational range of our FDI scheme as was defined in the
preceding case. It can also be observed that the percentage of the WMSNE is lower for the
cruise flight mode as compared to the climbing and landing modes due to the high variations of
the engine thrust.
E. Case 5: Concurrent Fault Detection, Isolation and Identification
In previous simulations only a single sensor fault is applied to the gas turbine engine. In this
case study, our proposed FDII scheme is used for concurrent fault scenarios that occur at different
stages of the flight profile. For this purpose, the hierarchical scheme described in Section IV-B
is used to diagnose the concurrent sensor faults. The fault detection times and the percentages
of the WMSNE as given by (23) are obtained for two selected scenarios with concurrent sensor
faults. The first scenario simulates the effects of a 6% bias fault applied to the TC sensor at
kf = 50 sec during the climbing mode and the concurrent 5% bias fault applied to the N sensor
at kf = 250 sec during the cruise mode. The estimated concurrent faults are depicted in Fig. 4.
The concurrent faults for the first scenario are detected and isolated after 1.7 and 2.4 seconds, respectively from their occurrence instants, whereas the percentages of the WMSNEs
are 0.0208% and 0.106%, respectively. The second scenario simulates the effects of a 4%
bias fault that is applied to the TT sensor at kf = 250 sec during the cruise mode and the
concurrent 6% bias fault that is applied to the PT sensor at kf = 450 sec during the landing
mode. The concurrent faults for the second scenario are detected and isolated after 0.9 and 0.4
seconds, respectively from their occurrence instants, whereas the percentages of the WMSNEs
are 0.0125% and 0.1702%, respectively.
F. Case 6: Comparison
In this section, the efficiency of our proposed MM-based sensor FDI scheme which utilizes
the MHKF as a detection filter is compared with different filtering methods including the linear
Kalman filter (LKF), extended Kalman filter (EKF), unscented Kalman filter (UKF), and the
cubature Kalman filter (CKF) [43]. It must be noted that in the following conducted experiments,
multiple combined residual vectors and covariance matrices associated with multiple operating
points are computed for both the HKF and LKF filtering methods to construct the multiple
HKF (MHKF) and multiple LKF (MLKF)-based FDI schemes, respectively. Moreover, our MMbased FDI scheme is compared with the interacting multiple model (IMM) approach [44], [45].
Different experiments are conducted in this section to compare the promptness of the fault
August 28, 2015

DRAFT

Sensor N Bias Fault

Sensor TC Bias Fault

22

40
20

Estimated Fault
Real Fault

0
-20
0

100

200

300

Time (sec)

400

500

600

(a)

600
400

Estimated Fault
Real Fault

200
0
-200
0

100

200

300

Time (sec)

400

500

600

(b)

Fig. 4: Estimated two concurrent sensor fault severities (a) a 6% bias fault that is applied to the
TC sensor at kf = 50 sec, and (b) a 5% bias fault that is applied to the N sensor at kf = 250
sec without considering any health parameter degradations.
detection and isolation scheme as well as the degree of robustness towards the health parameters
degradations among the above filtering methods. In all the experiments, similar process and
measurement noise signals are applied to the gas turbine engine and similar noise covariance
matrices are used to construct the above filters. The sampling rates for simulating the actual
engine and the nonlinear filtering approaches are 0.01 sec and 0.1 sec, respectively. Also, 2n + 1
sigma points, where n is the dimension of the state variable vector, are generated in the UKF
method with the tuning parameters α = 10−3 , β = 2 and κ = 0. For the CKF method, 2n
cubature points are √
generated through intersecting the unit sphere with the Cartesian axes and
are then scaled by n. Due to space limitations, the details regarding the formulations of the
above applied filtering methods are not included as they can be found in the literature.
Experiment 1: The 3% single bias fault is applied to different sensors at kf = 250 sec during
the cruise mode with no health parameter degradations. All the filters are capable of detecting
and isolating the injected fault correctly with different fault detection times. The sensor FDTs
are shown in TABLE XII for the detection filters that are utilized in the MM-based FDI scheme.
It can be concluded from TABLE XII that our proposed approach is capable of detecting and
isolating the sensor faults faster than other methods, particularly the MLKF method, given that the
OBEM can capture the nonlinear behaviour of the engine more effectively than others. Although,
August 28, 2015

DRAFT

23

1.2

1.2
P1

P2

P3

P4

P5

P

P6

Mode Probability

Mode Probability

0.5

0
0

100

200

300
400
Time (Sec)

500

P1

P2

P3

P4

P5

200

3

P

4

P

5

P

6

100

200

300

400

Time (Sec)

(d)

500

300
400
Time (Sec)

500

600

(c)

P

1

0.5

100

P

1.2

P6

Mode Probability

Mode Probability

2

(b)

1

0
0

P

0.5

0
0

600

(a)

1.2

1

1

1

600

1

P

2

P

3

P

4

P

5

P

6

0.5

0
0

100

200

300
400
Time(Sec)

500

600

(e)

Fig. 5: Mode probabilities for a 3% bias fault that is applied at kf = 250 sec to the TC sensor
by applying multiple (a) HKF, (b) LKF, (c) EKF, (d) UKF and (e) CKF based FDI approaches
with 2.5% degradations applied to the compressor health parameters of the actual engine.
TABLE XI: Absolute mean of residual signals associated with various types of nonlinear filters
corresponding to the healthy mode in the MM-based scheme for the healthy engine sensors
without health degradation.
Sensor
TC
PC
N
TT
PT

Absolute Mean of Residual Signal
MHKF
EKF
UKF
CKF
8.8329 × 10−4 0.00063 0.0058 0.0039
0.0044
0.0186 0.0127 0.0105
−4
4.0454 × 10
0.0032 0.0011 0.0009
8.1498 × 10−4
0.0081 0.0085 0.0069
0.0014
0.0028 0.0082 0.0074

the EKF, UKF, and CKF are well-known nonlinear filtering methods, however they detect and
isolate the injected sensor faults with more delay than our proposed method. According to our
observations that are included in TABLE XI, the HKF has smaller estimation error as compared
to the EKF, UKF, and CKF, which leads to faster FDI performance. It is also expected that
in case the OBEM is subjected to certain unmodeled dynamics with respect to a real engine,
the FDTs that are obtained with multiple EKF, UKF, and CKF based FDI schemes are more
comparable with our proposed MHKF-based FDI scheme.
Experiment 2: The robustness of our developed MHKF-based FDI approach is compared
with the other filtering methods here. The 3% bias fault is applied to the TC sensor at kf = 250
sec during the cruise mode with 2.5% degradations that are applied to the compressor health
parameters. In this experiment, it is assumed that the OBEM health parameters are not updated
August 28, 2015

DRAFT

24

and the 2.5% degradation is only applied to the compressor health parameters of the actual engine.
Therefore, the reference baselines associated with the compressor health parameters are one for
both the OBEM model and the nonlinear model that is used for deriving the other nonlinear
filtering methods. The mode probabilities that are generated in this scenario are depicted in Fig.
5. It can be concluded that our method outperforms the other approaches in terms of robustness
with respect to the health parameters degradations. It is also shown in Fig. 5 that the applied
fault can be detected by using the EKF and MLKF methods in addition to our method, but false
alarms are generated by these methods before the occurrence of the fault.
Experiment 3: Our proposed MHKF method is similar to the EKF in the sense of the
linearization process at certain operating points instead of the entire flight profile. Therefore,
one may expect that the EKF method can yield similar results to that of the MHKF method in
case the health parameters are also updated for the nonlinear engine model used for deriving
the EKF. In this experiment, the same engine model where its health parameters are updated is
used for the MM-based structures that utilize the MHKF and EKF schemes. The 3% bias fault is
applied to the TC sensor at kf = 250 sec during the cruise mode with λ = 0.99 and λ̂ = 0.975.
The estimated health parameters are updated similarly for both the EKF and MHKF schemes.
Fig. 6 shows the mode probabilities that are generated in this scenario. According to Fig. 6, the
EKF is also capable of correctly detecting and isolating the applied fault by using the health
parameters updating mechanism. However, the fault detection time and the computational time
of the EKF is far greater than that of the MHKF. The MM-based FDI method is capable of
detecting and isolating the applied fault in 2.7 sec and 20.8 sec by using the MHKF and EKF
methods, respectively. The same results are also obtained for the other fault scenarios, but are
not included here due to the space limitations.
Experiment 4: In this experiment, the performance of the FDI algorithm that utilizes our
proposed MHKF as well as the other filtering methods are compared in the interacting multiple
model (IMM) structure. For the IMM-based FDI structure, the sensor FDTs are computed for
the 3% single bias fault that is applied to the sensors with no health parameter degradations
as represented in the TABLE XIII. According to the simulation results, the fault is detected
sooner by using the IMM structure that utilizes all the detection filters except the MHKF. It
appears that the interactions among the multiple models do not improve the MHKF-based FDI
scheme performance. For simulating the IMM scheme which utilizes our proposed MHKF, the
off-diagonal elements of the transition probability matrix should be selected less than 0.001;
otherwise, the probability of the healthy mode monotonically decreases to reach almost 0.8
before the occurrence of a sensor fault which can enhance the risk of a false alarm. Therefore,
it is not recommended to use the IMM scheme for the MHKF method.
Experiment 5: In this experiment, the computational time of all the above filtering methods
are obtained and compared with our proposed MHKF approach in TABLE XIV. This comparison
confirms the advantages of implementing the MHKF approach in real-time applications. The
computational time is measured for only the filters corresponding to the healthy mode by the
use of “tic-toc” command in Matlab. It must be noted that the computational time of the MHKF
includes the time for not only simulating the healthy mode but also for running the OBEM
model. Therefore, based on the information provided in TABLE XIV, the MHKF scheme is
faster and more suitable than the other filtering methods to be utilized during the entire flight
August 28, 2015

DRAFT

25

TABLE XII: Sensor fault detection times using multiple non-interacting HKF, LKF, EKF, UKF
and CKF based FDI approaches without considering health parameter degradation.
MM-Based
FDI Approach
HKF
LKF
EKF
UKF
CKF

Detection Time (FDT) (sec)
PC
N
TT
PT
0.2
0.2
0.2
0.3
38.2 16.9 12.1
12.4
1
0.9
0.8
0.9
0.7
0.6
0.7
2.7
0.6
0.7
0.7
2.3

1.2

1.2

P

1

P

2

P

3

P

P

4

5

P

6

1

Mode Probability

Mode Probability

Fault
TC
0.2
8.3
1
0.6
0.6

0.5

0
0

100

200

300

400

Time (Sec)

500

600

(a)

1

P

1

P

2

P

3

P

4

P

5

P

6

0.5

0
0

100

200 300 400
Time (Sec)

500

600

(b)

Fig. 6: Mode probabilities for a 3% bias fault applied at kf = 250 sec to the TC sensor using
(a) the MHKF and (b) the EKF-based FDI methods, with λ = 0.99 and λ̂ = 0.975.

profile. Computing the Jacobian functions in each time step for the EKF method and also the
cubature and sigma points with their corresponding mean and covariance weight matrices for
the CKF and UKF methods, respectively are time consuming operations that lead to the higher
computational time compared to our proposed MHKF scheme. Moreover, the MLKF approach
computes the Kalman gain matrices on-line, whereas the MHKF utilizes the previously stored
gain matrices. However, if the same steady-state Kalman gain matrices are computed off-line

TABLE XIII: Sensor fault detection times using multiple interacting HKF, LKF, EKF, UKF and
CKF based FDI approaches without considering any health parameter degradation.
MM-Based
FDI Approach
HKF
LKF
EKF
UKF
CKF

August 28, 2015

Fault
TC
0.3
6.3
0.9
0.5
0.5

Detection Time (FDT) (sec)
PC
N
TT
PT
0.4
0.3
0.3
0.5
35.8 15.7 11.6
11.8
0.9
0.7
0.8
0.8
0.5
0.5
0.6
2.3
0.5
0.5
0.5
1.9

DRAFT

26

TABLE XIV: Computational time as measured in seconds for different filtering methods.
MHKF
131.21

MLKF
526.81

EKF
1184.2

UKF
2856.4

CKF
2673.44

and stored for use in the MLKF approach, the average computational time associated with the
healthy mode filters in the MLKF scheme is about 80 seconds, which is less than that of the
MHKF method. Under this situation, the extra computational time required by the MHKF is
associated with the OBEM operations. Therefore, the only disadvantage of the MHKF method
is the amount of memory that is needed for storing the look-up tables that include the matrices
Ai , C i and Kssi at different operating points.
Finally, in the following experiment results, it is shown that the number of required operating
points, and therefore the amount of required memory for the MHKF, is less than that of the
MLKF-based FDI method.
Experiment 6: In this experiment, the number of required operating points is investigated
by measuring the mean of the combined residual signals in presence of the health parameter
degradations for the MLKF and MHKF methods. To cover the entire flight profile, five operating
points are generated for the MHKF method as described in TABLE II. Note that the MLKF
method still needs more operating points since it is not capable of tracking rapid variations of
the thrust during the climbing and landing modes, and therefore it generates false alarms due to
the higher estimation error as compared to the MHKF scheme. TABLE XV compares the mean
of the MHKF and MLKF combined residual signals that are associated with the healthy mode
as a function of the number of operating points. Fig. 7 displays the mode probabilities that are
generated for the MHKF and MLKF methods for different number of operating points. In this
experiment, the engine has no sensor fault, but a 3% degradation is applied to the actual engine
compressor health parameters. As stated in Section VI-A, the MHKF is robust with respect to the
applied health degradations even without updating the OBEM model health parameters, although
the MHKF needs to be updated for larger degradation levels. Based on the results obtained, it
can be concluded that the mean of the combined residual signal associated with the healthy mode
can be reduced by increasing the number of the operating points. It is expected that the MLKF
estimation error decreases further by increasing the number of the operating points. However,
we are not capable of perfectly avoiding generation of a false alarm by using the MLKF method
corresponding to our particular flight profile by increasing the number of the operating points,
since the MLKF cannot track high variations of the thrust in our application. Therefore, it may
be feasible to apply the MLKF with the combined PWL models if the engine thrust increases
with the slower speed than ours since the PWL model weights calculation is a time-consuming
process which may diminish the ability of the MLKF to track the rapid variations of the thrust.
Therefore, it can be concluded that the MM-based FDI scheme which utilizes our proposed
MHKF approach is capable of promptly detecting and isolating various sensor bias faults during
the entire flight profile by updating the OBEM health parameters reference baselines and even
without updating the OBEM model for certain levels of health parameter degradations. Our
scheme also estimates the engine state and output variables with more accuracy and smaller
number of operating points as compared to that by the MLKF-based scheme.

August 28, 2015

DRAFT

27

1.2

P1

P2

P3

P4

P5

1.2

P6

P

1

2

P

3

P

4

P

5

P

1.2

100

200

300
400
Time (Sec)

500

600

1

P

2

P

3

P

4

P

5

P

6

1

0.5

0
0

P

6

Mode Probability

0.5

0
0

P

1
Mode Probability

Mode Probability

1

100

200

(a)

300
400
Time (Sec)

500

0.5

0
0

600

100

200

300
400
Time (Sec)

500

600

500

600

(c)

(b)
1.2
P

1

P

2

P

3

P

4

P

5

P

6

Mode Probability

1

0.5

0
0

(d)

(e)

100

200

300
400
Time (Sec)

(f)

Fig. 7: Mode probabilities that are generated for the MM-based FDI method by utilizing (a-c)
the MHKF and (d-f) the MLKF with one, three and five operating points.
TABLE XV: Mean of the combined residual signals (γ 1c ) associated with the healthy mode for
both the MHKF and MLKF methods in terms of the number of the operating points, whereas
3% degradation is applied to the compressor health parameters.
Sensor
TC
PC
N
TT
PT

1 operating point
HKF
LKF
0.109
0.2305
0.2286
0.943
0.0255 0.0123
0.0813 0.3339
0.0152 0.0477

3 operating points
MHKF MLKF
0.077
0.1091
0.0771
0.1523
0.0167
0.0053
0.0745
0.2372
0.0116
0.036

5 operating points
MHKF MLKF
0.0119
0.0857
0.0216
0.1315
0.0003
0.0014
0.0197
0.1008
0.0006
0.0049

VII. CONCLUSIONS
In this paper, a novel approach for both single and concurrent sensor fault detection, isolation
and identification/estimation (FDII) for gas turbine engines is proposed. Our methodology is
based on the developed Hybrid Kalman Filter (HKF) as the detection filter of a hierarchical
multiple-model based structure. Despite the use of linear Kalman filters, the HKF is capable
of capturing the nonlinearities of the system by integrating a nonlinear on-board engine model
(OBEM) with piecewise linear (PWL) models to cover the entire operating range of the engine.
Compared to the multiple linear Kalman filter (MLKF), our proposed approach requires fewer
number of operating points although each corresponds to a larger operating range.
Another important contribution of this work is in the inclusion of the effects of health parameter
degradations in our proposed sensor FDII scheme through updating the OBEM health parameters
reference baselines that enables one to prevent false alarms and incorrect fault detections.
August 28, 2015

DRAFT

28

Therefore, unlike most of the previous work in the literature that have not considered the
influence of health parameters variations on the performance of the developed FDII approaches,
our proposed sensor FDII scheme is sufficiently more effective for use in the entire flight profile.
It is furthermore suitable during the engine life cycle by updating the OBEM health parameters
reference baselines, and by replacing the steady-state values with the OBEM states and outputs
and integration of the designed PWL models. Moreover, our proposed MM-based sensor FDI
scheme is integrated with the modified GLR method to estimate a sensor fault severity. Finally,
it is demonstrated through extensive simulation studies that the accuracy and robustness of
our proposed MHKF-based FDI approach against the engine health parameter degradations are
significantly superior to that of those generated by other investigated methods such as the MLKF,
EKF, UKF, and CKF.
A PPENDIX A
In Section II, for sake of notational simplicity the same notations were used for both versions
of the OBEM. However, here we use XNOBEM and YNOBEM to designate the nominal OBEM
(λ = 1) that is equivalent to (2), and XOBEM and YOBEM to designate the OBEM with the
updated health parameters (λ̂ 6= 1) that is equivalent to (5).
Let the nominal continuous-time OBEM model (NOBEM) be represented by:
ẊNOBEM (t) = fc (XNOBEM (t), HNOBEM (t), U (t)),
YNOBEM (t) = gc (XNOBEM (t), HNOBEM (t)).
To construct the HKF, the nominal continuous-time OBEM is initially linearized at a certain
operating point (Xss , Uss , Yss ) without considering degradations as follows:
∂fc
∂fc
|U (U (t) − Uss ) + ∆fc ,
|Xss (XNOBEM (t) − Xss ) +
∂XNOBEM
∂U ss
∂gc
YNOBEM (t) = gc (Xss , Uss ) +
|X (XNOBEM (t) − Xss ) + ∆gc ,
∂XNOBEM ss

ẊNOBEM (t) = fc (Xss , Uss ) +

where fc (Xss , Uss ) ≡ 0 and gc (Xss , Uss ) ≡ Yss . The health parameters are the state-dependent
functions whose effects are incorporated in the above partial derivatives. The linear continuoustime model is now discretized by using a sufficiently small sampling period to yield the corresponding discrete-time linear model:
XNOBEM (k + 1) − Xss = A|Xss (XNOBEM (k) − Xss ) + B|Uss (U (k) − Uss ) + ∆F,
YNOBEM (k) = Yss + C|Xss (XNOBEM (k) − Xss ) + ∆G,
RT
where A|Xss = exp(Ac |Xss Ts ), B|Uss = ( 0 s exp(Ac |Xss τ )dτ )Bc |Uss , and C|Xss = Cc |Xss are
obtained in terms of the state-space matrices that are associated with the continuous-time linear
∂gc
∂fc
c
|Xss , Bc |Uss = ∂f
| , and Cc |Xss = ∂XNOBEM
|Xss , and the sampling
model Ac |Xss = ∂XNOBEM
∂U Uss
period is denoted by Ts . Also, ∆fc , ∆gc , ∆F and ∆G represent the higher order terms in
the linearization process. To construct the HKF, it is necessary to update the OBEM health
parameters. Therefore, the nominal model is only used once to derive the state-space matrices
and the steady-state values.
Assumption 1. The linearization of the nominal OBEM with λ = 1 and the updated OBEM
with λ̂ 6= 1 generate approximately the same matrices and steady-state values.
August 28, 2015

DRAFT

29

Therefore, given that Assumption 1 holds, the linearized version of the updated OBEM can
be derived as follows:
XOBEM (k + 1) − Xss ≈ A|Xss (XOBEM (k) − Xss ) + B|Uss (U (k) − Uss ) + ∆F,
YOBEM (k) ≈ Yss + C|Xss (XOBEM (k) − Xss ) + ∆G.

(24)

Moreover, the relationship between the instantaneous values of XOBEM (k) and YOBEM (k) and
the steady-state values associated with each linearized model can be given by:
XOBEM (k) = Xss + Xl (k),
YOBEM (k) = Yss + Yl (k),

(25)

where Xl (k) and Yl (k) represent the perturbations from the steady-state values as well as
variations of the OBEM state variables and outputs due to updating the health parameters
reference baselines. The linear Kalman filter can be designed as follows:
X̂(k + 1) − Xss = A|Xss (X̂(k) − Xss ) + B|Uss (U (k) − Uss ) + Kss (Y (k) − Ŷ (k)),
Ŷ (k) = Yss + C|Xss (X̂(k) − Xss ),
which can be rewritten by using (25) as follows:
X̂(k + 1) − XOBEM (k) + Xl (k) = A|Xss (X̂(k) − XOBEM (k)) + Kss (Y (k) − Ŷ (k)) + A|Xss Xl (k)
+ B|Uss (U (k) − Uss ),
Ŷ (k) = YOBEM (k) + C|Xss (X̂(k) − XOBEM (k)) + C|Xss Xl (k) − Yl (k).

(26)

Using (24), the last two terms in the R.H.S of (26) can be rewritten as follows:
A|Xss (XOBEM (k) − Xss ) + B|Uss (U (k) − Uss ) ≈ XOBEM (k + 1) − Xss − ∆F,
C|Xss (XOBEM (k) − Xss ) − Yl (k) ≈ YOBEM (k) − Yss − Yl (k) − ∆G = −∆G.
Therefore, the discrete-time form of the HKF can be obtained as follows:
X̂(k + 1) − XOBEM (k + 1) ≈ A|Xss (X̂(k) − XOBEM (k)) + Kss (Y (k) − Ŷ (k)) − ∆F,
Ŷ (k) ≈ YOBEM (k) + C|Xss (X̂(k) − XOBEM (k)) − ∆G,
which will lead to (6) for the particular operating point if the higher order terms ∆F and ∆G
as well as the approximation error due to Assumption 1 are neglected.
R EFERENCES
[1] R. Isermann and P. Balle, “Trends in the application of model-based fault detection and diagnosis of technical processes,”
Control Engineering Practice, vol. 5, no. 5, pp. 709–719, 1997.
[2] P. M. Frank, “Fault diagnosis in dynamic systems using analytical and knowledge-based redundancy: A survey and some
new results,” Automatica, vol. 26, no. 3, pp. 459–474, 1990.
[3] W.-M. Lin and T.-C. Ou, “Unbalanced distribution network fault analysis with hybrid compensation,” Generation,
Transmission & Distribution, IET, vol. 5, no. 1, pp. 92–100, 2011.
[4] X. Zhang, L. Tang, and J. Decastro, “Robust fault diagnosis of aircraft engines: a nonlinear adaptive estimation-based
approach,” IEEE Trans. on Control Systems Technology, vol. 21, no. 3, pp. 861–868, 2013.
[5] V. Venkatasubramanian, R. Rengaswamy, K. Yin, and S. N. Kavuri, “A review of process fault detection and diagnosis:
Part I: Quantitative model-based methods,” Computers & Chemical Engineering, vol. 27, no. 3, pp. 293–311, 2003.

August 28, 2015

DRAFT

30

[6] I. Hwang, S. Kim, Y. Kim, and C. E. Seah, “A survey of fault detection, isolation, and reconfiguration methods,” IEEE
Trans. on Control Systems Technology, vol. 18, no. 3, pp. 636–653, 2010.
[7] J. C. d. Silva, A. Saxena, E. Balaban, and K. Goebel, “A knowledge-based system approach for sensor fault modeling,
detection and mitigation,” Expert Systems with Applications, vol. 39, no. 12, pp. 10977–10989, 2012.
[8] T. E. Menke and P. S. Maybeck, “Sensor/actuator failure detection in the vista F-16 by multiple model adaptive estimation,”
IEEE Trans. on Aerospace and Elec. Sys., vol. 31, no. 4, pp. 1218–1229, 1995.
[9] P. S. Maybeck, “Multiple model adaptive algorithms for detecting and compensating sensor and actuator/surface failures
in aircraft flight control systems,” Int. J. of Robust and Nonlinear Control, vol. 9, no. 14, pp. 1051–1070, 1999.
[10] N. Meskin, E. Naderi, and K. Khorasani, “Fault diagnosis of jet engines by using a multiple model-based approach,” ASME
Turbo Expo, no. GT2010-23442, pp. 319–329, 2010.
[11] N. Meskin, E. Naderi, and K. Khorasani, “A multiple model-based approach for fault diagnosis of jet engines,” IEEE
Trans. on Control Sys. Tech., pp. 254–262, 2013.
[12] E. Naderi, N. Meskin, and K. Khorasani, “Nonlinear fault diagnosis of jet engines by using a multiple model-based
approach,” J. of Eng. for Gas Turbines and Power, vol. 134, no. 1, p. 011602, 2012.
[13] B. Pourbabaee, N. Meskin, and K. Khorasani, “Sensor fault detection and isolation using multiple robust filters for linear
systems with time-varying parameter uncertainty and error variance constraints,” in IEEE Conf. on Control Applications
(CCA), pp. 382–389, 2014.
[14] A. J. Volponi, “Foundation of gas path analysis (part I and II),” Gas Turbine Condition Monitoring and Fault Diagnosis,
vol. 1, 2003.
[15] T. Kobayashi and D. L. Simon, “Evaluation of an enhanced bank of Kalman filters for in-flight aircraft engine sensor fault
diagnostics,” J. of Eng. for Gas Turbines and Power, vol. 127, p. 497, 2005.
[16] T. Kobayashi and D. L. Simon, “Application of a bank of Kalman filters for aircraft engine fault diagnostics,” in ASME
Turbo Expo, collocated with the International Joint Power Generation Conference, pp. 461–470, 2003.
[17] T. Kobayashi and D. L. Simon, “Hybrid Kalman filter approach for aircraft engine in-flight diagnosis: Sensor fault detection
case,” J. of Eng. for Gas Turbines and Power, vol. 129, no. 3, pp. 746–754, 2007.
[18] S. Simani, C. Fantuzzi, and S. Beghelli, “Diagnosis techniques for sensor faults of industrial processes,” IEEE Trans. on
Control Systems Technology, vol. 8, no. 5, pp. 848–855, 2000.
[19] A. J. Volponi, H. DePold, R. Ganguli, and C. Daguang, “The use of Kalman filter and neural network methodologies in
gas turbine performance diagnostics: a comparative study,” J. of Eng. for Gas Turbines and Power, vol. 125, pp. 917–924,
2003.
[20] R. Ganguli, “Application of fuzzy logic for fault isolation of jet engines,” J. of Eng. for Gas Turbines and Power, vol. 125,
pp. 617–623, 2003.
[21] S. Sampath and R. Singh, “An integrated fault diagnostics model using genetic algorithm and neural networks,” J. of Eng.
for Gas Turbines and Power, vol. 128, pp. 49–56, 2006.
[22] S. Rahme and N. Meskin, “Adaptive sliding mode observer for sensor fault diagnosis of an industrial gas turbine,” Control
Engineering Practice, vol. 38, pp. 57–74, 2015.
[23] E. Tsoutsanis, N. Meskin, M. Benammar, and K. Khorasani, “Transient gas turbine performance diagnostics through
nonlinear adaptation of compressor and turbine maps,” J. of Engineering for Gas Turbines and Power, vol. 137, no. 9,
p. 091201, 2015.
[24] E. Tsoutsanis, N. Meskin, M. Benammar, and K. Khorasani, “A component map tuning method for performance prediction
and diagnostics of gas turbine compressors,” Applied Energy, vol. 135, pp. 572–585, 2014.
[25] L. Shang and G. Liu, “Sensor and actuator fault detection and isolation for a high performance aircraft engine bleed air
temperature control system,” IEEE Trans. on Control Systems Technology, vol. 19, no. 5, pp. 1260–1268, 2011.
[26] R. Babuska and H. B. Verbruggen, “Fuzzy set methods for local modeling and identification,” 1997.
[27] X. C. Zhang, A. Visala, A. Halme, and P. Linko, “Functional state modelling approach for bioprocesses: local models for
aerobic yeast growth processes,” J. of Process Control, vol. 4, no. 3, pp. 127–134, 1994.
[28] M. Meila and M. I. Jordan, “Markov mixtures of experts,” Multiple Model Approaches to Modeling and Control, pp. 145–
166, 1997.
[29] M. M. Seron, J. A. De Dona, and J. H. Richter, “Integrated sensor and actuator fault-tolerant control,” Int. J. of Control,
vol. 86, no. 4, pp. 689–708, 2013.
[30] A. S. Willsky and H. L. Jones, “A generalized likelihood ratio approach to the detection and estimation of jumps in linear
systems,” IEEE Trans. on Automatic Control, vol. 21, no. 1, pp. 108–112, 1976.
[31] J. Prakash, S. Narasimhan, and S. C. Patwardhan, “Integrating model based fault diagnosis with model predictive control,”
Industrial & Engineering Chemistry Research, vol. 44, no. 12, pp. 4344–4360, 2005.
[32] S. C. Patwardhan, S. Manuja, S. Narasimhan, and S. L. Shah, “From data to diagnosis and control using generalized
orthonormal basis filters. part II: Model predictive and fault tolerant control,” J. of Process Control, vol. 16, no. 2,
pp. 157–175, 2006.
[33] P. M. Frank and X. Ding, “Survey of robust residual generation and evaluation methods in observer-based fault detection
systems,” J. of process control, vol. 7, no. 6, pp. 403–424, 1997.

August 28, 2015

DRAFT

31

[34] R. Seliger and P. M. Frank, “Fault-diagnosis by disturbance decoupled nonlinear observers,” in Proceedings of the 30th
IEEE Conf. on Decision and Control, pp. 2248–2253, 1991.
[35] D. Simon, “A comparison of filtering approaches for aircraft engine health estimation,” J. of Aerospace Science and
Technology, vol. 12, no. 4, pp. 276–284, 2008.
[36] B. Pourbabaee, N. Meskin, and K. Khorasani, “Multiple-model based sensor fault diagnosis using hybrid Kalman filter
approach for nonlinear gas turbine engines,” in American Control Conference (ACC), pp. 4717–4723, 2013.
[37] R. Murray-Smith and T. Johansen, Multiple Model Approaches to Nonlinear Modelling and Control. CRC press, 1997.
[38] X. R. Li and Y. Zhang, “Numerically robust implementation of multiple-model algorithms,” IEEE Trans. on Aerospace
and Electronic Systems, vol. 36, no. 1, pp. 266–278, 2000.
[39] J. Prakash, S. C. Patwardhan, and S. Narasimhan, “A supervisory approach to fault-tolerant control of linear multivariable
systems,” Industrial & Engineering Chemistry Research, vol. 41, no. 9, pp. 2270–2281, 2002.
[40] A. P. Deshpande and S. C. Patwardhan, “Online fault diagnosis in nonlinear systems using the multiple operating approach,”
Ind. Eng. Chem. Res., vol. 47, no. 17, pp. 6711–6726, 2008.
[41] N. L. R. Group, “GSP software,” http://www.gspteam.com, 2011.
[42] T. Fawcett, “An introduction to ROC analysis,” Pattern recognition letters, vol. 27, no. 8, pp. 861–874, 2006.
[43] I. Arasaratnam and S. Haykin, “Cubature Kalman filters,” IEEE Trans. on Automatic Control, vol. 54, no. 6, pp. 1254–1269,
2009.
[44] Y. Zhan and J. Jiang, “An interacting multiple-model based fault detection, diagnosis and fault-tolerant control approach,”
in Proceedings of the 38th IEEE Conf. on Decision and Control, vol. 4, pp. 3593–3598, 1999.
[45] J. Ru and X. R. Li, “Variable-structure multiple-model approach to fault detection, identification, and estimation,” IEEE
Trans. on Control Systems Technology, vol. 16, no. 5, pp. 1029–1038, 2008.

August 28, 2015

DRAFT

