- Email: [email protected]
Sufficiently Exciting Inputs for Structurally Identifiable Systems Biology Models⁎
Proceedings, 7th IFAC Conference on Proceedings, 7th IFAC Conference on Foundations of Systems Biology inon Engineering Proceedings, 7th IFAC Conference Available online at www.sciencedirect.com Foundations of Systems Biology inon Engineering Proceedings, 7th IFAC August Conference Chicago, Illinois, USA, 5-8,in2018 Foundations of Systems Biology Engineering Chicago, Illinois, USA, August 5-8, 2018 Foundations of Systems Biology in Engineering Chicago, Illinois, USA, August 5-8, 2018 Chicago, Illinois, USA, August 5-8, 2018
ScienceDirect
IFAC PapersOnLine 51-19 (2018) 16–19
Sufficiently Exciting Inputs for Structurally Sufficiently Exciting Inputs for Structurally Sufficiently Exciting Inputs for Structurally Identifiable Systems Biology Models Sufficiently Exciting Inputs for Structurally Identifiable Systems Biology Identifiable Systems Biology Models Models Identifiable Systems Biology Models ∗ Alejandro F. Villaverde Neil D. Evans ∗∗
∗ ∗∗ Alejandro F. Villaverde ∗∗ ∗ Neil D. Evans∗∗∗ Michael J. Julio Alejandro F.Chappell Villaverde NeilR. D.Banga Evans∗∗∗ ∗∗ ∗ Michael J. Chappell Julio R. Banga Alejandro Villaverde NeilR. D.Banga Evans∗ ∗∗ Julio Michael J.F.Chappell ∗∗ ∗ Michael J. Chappell Julio R. Banga ∗ ∗ Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain ∗ Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain (e-mail: ∗ (e-mail: [email protected]). [email protected]). ∗∗ Bioprocess Engineering Group, IIM-CSIC, 36208,CV4 Spain University of Warwick,Vigo Coventry 7AL, (e-mail: [email protected]). ∗∗ School of Engineering, University of Warwick, Coventry CV4 7AL, (e-mail: [email protected]). ∗∗ School of Engineering, UK School of Engineering, University of Warwick, Coventry CV4 ∗∗ UK of Warwick, Coventry CV4 7AL, School of Engineering, University 7AL, UK UK Abstract: its value can theoretically be estimated Abstract: A A parameter parameter is is structurally structurally identifiable identifiable if if its value can theoretically be estimated by observing the model output. Structural identifiability is a desirable property be in estimated biological Abstract: A parameter is structurally identifiable if its value can theoretically by observingA the model output. Structural identifiability is aa desirable property in biological Abstract: parameter is structurally identifiable if its value can theoretically be by observingif the model output. Structural identifiability is its desirable property in estimated biological modelling: a parameter is structurally unidentifiable, estimated numerical value modelling: if a parameter is structurally unidentifiable, its estimated numerical value is is by observingifand the model predictions output. Structural identifiability is its a desirable property in biological modelling: a parameter is structurally unidentifiable, estimated numerical value is meaningless, model of unmeasured state variables can be wrong, compromising meaningless, and model predictions of unmeasured state variables can be wrong, compromising modelling: if a parameter is structurally unidentifiable, its estimated numerical value is the ability the model to insight. Structural on the meaningless, of unmeasured state variablesidentifiability can be wrong,depends compromising the ability of of and the model model predictions to provide provide biological biological insight. Structural identifiability depends on the meaningless, and model predictions of unmeasured state variables can be wrong, compromising system dynamics, observation function (model output), initial and the ability of the model to provide biological Structural identifiability depends inputs. on the system dynamics, observation function (modelinsight. output), initial conditions, conditions, and external external inputs. the ability of we thefocus model to provide biological insight. Structural identifiability depends on the system dynamics, observation function (model output), initial conditions, andanalysis external inputs. In this paper on the last factor. Methods for structural identifiability typically In this paper we focus on the last factor.(model Methods for structural identifiability analysis typically system dynamics, observation function output), initial conditions, and external inputs. classify a model as identifiable provided that it is fed with sufficiently exciting inputs. For In this paper we focus on the last factor. Methods for structural identifiability analysis typically classify a model as identifiable provided that it for is structural fed with sufficiently For exciting inputs. In this paper we focus onmay the last factor. Methods identifiability analysis typically example, a given model require a time-varying input to be structurally identifiable, while classify a model as identifiable provided that it is fed with sufficiently exciting inputs. For example,a amodel given model may require a time-varying input to besufficiently structurally identifiable, while classify as identifiable provided that it is fed with exciting inputs. For for another model a constant non-zero input may be enough. Here we present a method that example, a given model may require a time-varying input to be structurally identifiable, while for another modelmodel a constant non-zero input may beinput enough. Here we present aa method that example, a given may require a time-varying to be structurally identifiable, while for another model a constant non-zero input may be enough. Here we present method that determines how sufficiently exciting an should be. The builds on STRIKEdetermines how sufficiently exciting an input input The approach approach onathe the STRIKEfor another model a constant non-zero input should may bebe. enough. Here we builds present method that GOLDD toolbox, which considers structural identifiability as generalized observability. The determines how sufficiently exciting an input should be. The approach builds on the STRIKEGOLDD toolbox, which considers structural identifiability as generalized observability. The determines how sufficiently exciting an input should be. The approach builds on the STRIKEapproach incorporates extended Lie derivatives, which correctly assess structural identifiability GOLDD toolbox, which considers structural identifiability as generalized observability. The approach incorporates extended Lie structural derivatives,identifiability which correctly assess structural identifiability GOLDD toolbox, which considers as generalized observability. The in the case of time-varying inputs. The procedure can also be used to determine the type approach incorporates extended Lie derivatives, which correctly assess structural identifiability in the case of time-varying inputs. The procedure cancorrectly also be used tostructural determineidentifiability the type of of approach incorporates extended Lie derivatives, which assess in the case of time-varying inputs. The procedure can also be used to determine the type of input profile that is to the identifiable. This capability is helpful when input profile is required required inputs. to make makeThe the parameters parameters identifiable. that This capability is helpful when in the case of time-varying procedure can also be used to determine the type of designing new experiments for the purpose of parameter estimation. input profile that is required to make the parameters identifiable. This capability is helpful when designing new experiments for the purpose of parameter estimation. input profile that is required to make the parameters identifiable. This capability is helpful when designing new experiments for the purpose of parameter estimation. © 2018, IFAC of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. designing new(International experimentsFederation for the purpose of parameter estimation. Keywords: properties, Observability, Keywords: Identifiability, Identifiability, Structural Structural Identifiability, Identifiability, Structural Structural properties, Observability, Keywords: Identifiability, Structural Identifiability, Structural properties, Observability, Parameter identification, estimation, Parameter identification, Parameter Parameter estimation, Inputs. Inputs. Keywords: Identifiability, Structural Identifiability, Structural properties, Observability, Parameter identification, Parameter estimation, Inputs. Parameter identification, Parameter estimation, Inputs. 1. INTRODUCTION Molenaar, 2015). There are a number of software tools 1. tools Molenaar, There are aa number software 1. INTRODUCTION INTRODUCTION Molenaar, 2015). 2015). There are number of of methodologies, software tools implementing some of the aforementioned implementing some of the aforementioned methodologies, 1. INTRODUCTION Molenaar, 2015). There are a number of software tools A parameter is structurally identifiable if it is theoreti- implementing such as DAISY (Bellu et al., 2007), GenSSI (Chi¸ ss et al., some of the aforementioned methodologies, A parameter is structurally identifiable if it is theoretisuch as DAISY (Bellu et al., 2007), GenSSI (Chi¸ et al., implementing some of the aforementioned methodologies, A parameter is structurally identifiable if it is theoretically possible to determine its true value from noiseless such 2011), EAR (Karlsson 2012), COMBOS (Meshkat as DAISY (Bellu et al., 2007), GenSSI (Chi¸ s et al., cally possible to determine its true value from noiseless 2011), EAR (Karlsson et al., 2012), COMBOS (Meshkat A parameter is structurally identifiable if it is theoretias2014), DAISY (Bellu et al., 2007), GenSSI (Chi¸ s 2016), et al., cally possible and to determine true value noiseless such data (Cobelli DiStefano,its1980). Since from identifiability 2011), EAR (Karlsson et al., 2012), COMBOS (Meshkat et al., STRIKE-GOLDD (Villaverde et al., data and DiStefano, 1980). Since identifiability et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), cally (Cobelli possible to other determine true noiseless EAR (Karlsson et al.,et2012), COMBOS (Meshkat data (Cobelli and DiStefano, 1980). Since from identifiability depends, among things,itson the value system’s inputs and 2011), and Data2Dynamics (Raue al., 2015). These methods et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), depends, among other things, on the system’s inputs and and Data2Dynamics (Raue et al., 2015). These methods data (Cobelli and DiStefano, 1980). Since identifiability et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), depends, among other things, on the system’s inputs and outputs, it can vary as a result of the experiment design. and assess whether a given input-output configuration is in Data2Dynamics (Raue et al., 2015). These methods outputs, it can as a the experiment design. whether a given input-output configuration is in depends, other things, the and assess and Data2Dynamics (Raue et al., 2015). These methods outputs, can vary vary as a result resultonof of thesystem’s experiment The termitamong qualitative experiment design refers toinputs thedesign. selecprinciple sufficient to obtain a structurally identifiable assess whether a given input-output configuration is in The term qualitative experiment design refers to the selecprinciple sufficient to obtain a structurally identifiable outputs, it can vary as a result of the experiment design. assess whether a given input-output configuration is in The qualitative experiment design to the selection term of input and output ports in orderrefers to maximize the principle sufficient to Cobelli, obtain a1992). structurally identifiable model (Saccomani and Algorithms that tion of input and output ports in order to maximize the model (Saccomani and Cobelli, Algorithms that auauThe term experiment design to the selecsufficient tothe obtain a1992). structurally identifiable tion of input and output ports in orderrefers to maximize the principle number ofqualitative structurally identifiable parameters (Walter and tomate the search for subsets of outputs that guarantee model (Saccomani and Cobelli, 1992). Algorithms that aunumber of structurally identifiable parameters (Walter and tomate the search for the subsets of outputs that guarantee tion of input and output ports in order to maximize the model (Saccomani and Cobelli, 1992). Algorithms that aunumber of structurally identifiable parameters (Walter and Pronzato, 1990). structural identifiability are also available (August, 2009; tomate the search for the subsets of outputs that guarantee Pronzato, identifiability are also ofavailable (August, 2009; number of 1990). structurally identifiable parameters (Walter and structural tomate the search for the subsets outputs that guarantee Pronzato, 1990). Anguelova et identifiability Some structural identifiability methods, such as the simi- structural Anguelova et al., al., 2012). 2012). are also available (August, 2009; Pronzato, 1990). identifiability structural Some structural methods, such as the simiAnguelova identifiability et al., 2012). are also available (August, 2009; Some structural identifiability methods, such as the similarity transformation approach (Evans et al., 2002; Yates Anguelova Determining etthe al.,necessary 2012). set larity transformation approach (Evans et al., 2002; the necessary set of of inputs inputs and and outputs outputs is is only only Some identifiability the Yates simi- Determining larity transformation approach (Evans etsuch al., as 2002; Yates et al.,structural 2009) and direct testmethods, (Denis-Vidal and Jolypart of experiment design. Additionally, it is necessary Determining the necessary set of inputs and outputs is only et al., 2009) and direct test (Denis-Vidal and Jolypart of experiment design. Additionally, it is necessary larity transformation approach (Evans et al., 2002; Yates Determining the necessary setAdditionally, ofnumber inputs and it outputs is only et al., 2009) andaredirect test (Denis-Vidal Joly- part Blanchard, 2000) applicable to autonomousand systems, of experiment design. is necessary to aspects as of Blanchard, 2000) aredirect applicable to systems, to define define aspects such such as the theAdditionally, number and andittiming timing of the the et systems al., 2009) and test approaches, (Denis-Vidal and Joly- part of experiment design. is necessary Blanchard, 2000) applicable to autonomous autonomous systems, i.e. with noare input. Other such as those samples, the system’s initial state, and characterization of to define aspects such as the number and timing of the i.e. systems with no input. Other approaches, such as those samples, the system’s initial state, and characterization of Blanchard, 2000) are applicable to autonomous systems, to define aspects such as the number and timing of the i.e. systems with no input. Other approaches, such as those based on power series (Pohjanpalo, 1978), differential al- samples, the inputs. This task is known as quantitative experiment the system’s initial state, and characterization of based on series (Pohjanpalo, 1978), differential althe inputs. This task isinitial known as quantitative experiment i.e. systems with noGlad, input. Other approaches, such as those samples, the system’s state, and characterization of based on power power (Pohjanpalo, 1978), differential al- design, gebra (Ljung andseries 1994; Bellu et al., 2007), implicit and aims at maximizing practical identifiability. the inputs. This task is known as quantitative experiment gebra andseries Glad,(Pohjanpalo, 1994; Bellu et al., implicit andThis aimstask at maximizing practical identifiability. based (Ljung on power differential al- design, the the inputs. is we known as quantitative experiment gebra (Ljung 1994; Bellu et al., 2007), 2007), implicit functions (Xiaand andGlad, Moog, 2003), or 1978), differential geometry design, and aims at maximizing practical identifiability. In present paper address a different, although functions (Xia and Moog, 2003), or differential geometry In the present we addresspractical a different, although gebra (Ljung 1994; Bellu al.,to2007), implicit and aimspaper at maximizing identifiability. functions (Xia andGlad, Moog, 2003), or et differential geometry (Villaverde etand al., 2016), can be applied systems with design, related, problem. We seek to determine, a priori, what In the present paper we address a different, although (Villaverde et al., 2016), can be applied to systems with related, problem. We seek to determine, a priori, what functions (Xia and Moog, 2003), or differential geometry In the present paper we address a different, although (Villaverde et al., 2016), can be applied to systems with external inputs. Additionally, it is possible (for a given set related, mathematical form the inputs must have in order to guarproblem. We seek to determine, a priori, what external inputs. Additionally, it is possible (for a given set mathematical form the inputs must have in order to guar(Villaverde et al., 2016), can be applied to systems with problem. We seek to must determine, aorder priori, what external inputs.with Additionally, it is possible a given set related, of experiments defined inputs) to assess(for identifiability mathematical form the inputs have in to guarantee structural identifiability. For example, for a given of experiments defined to assess(for identifiability structural identifiability. For example, for a given external inputs.with Additionally, it is possible a given set antee mathematical form the inputsconstant must have in order to guarof experiments with defined inputs) inputs) to identifiability with numerical approaches based onassess profile likelihoods antee structural identifiability. For example, for a given model structure a non-zero input may not be with numerical approaches based on profile likelihoods model structural structure identifiability. a non-zero constant input may be of experiments with defined inputs) assess identifiability For example, for anot given with numerical approaches based to onmatrix profile likelihoods (Raue et al., 2009) or the sensitivity (Stigter and antee sufficiently exciting and a ramp input may be necessary. To model structure a non-zero constant input may not be (Raue et al., 2009) or the sensitivity and sufficiently exciting and a ramp input may be necessary. To with numerical approaches based onmatrix profile(Stigter likelihoods model structure a non-zero constant input may not be (Raue et al., 2009) or the sensitivity matrix (Stigter and aa certain extent, this situation resembles the relationship This project has received funding from the European Union’s Horisufficiently exciting and a ramp input may be necessary. To certain extent, this situation resembles thenecessary. relationship (Raue et al.,has 2009) or the sensitivity matrix (Stigter and sufficiently exciting and a ramp input may be To This project received funding from the European Union’s Horia certain extent, this situation resembles the relationship between structural identifiability and conditions: zon 2020 research innovation programme under grant agreement This project hasand received funding from the European Union’s Horibetween identifiability and initial initial conditions: zon 2020 research innovation programme under grant agreement a certain structural extent, this situation resembles theidentifiability relationship This project hasand received funding from the European Union’s Horigenerally, methods that analyse structural between structural identifiability and initial conditions: No 686282 (“CANPATHPRO”). zon 2020 research and innovation programme under grant agreement generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). between structural identifiability and initial conditions: zon 2020 research and innovation programme under grant agreement generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). 2405-8963 © 2018, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved.
Copyright © 2018 IFAC 16 Peer review©under responsibility of International Federation of Automatic Copyright 2018 IFAC 16 Control. Copyright © 2018 IFAC 16 10.1016/j.ifacol.2018.09.015 Copyright © 2018 IFAC 16
IFAC FOSBE 2018 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
yield results that are valid for almost all initial conditions; however, a model classified as structurally identifiable may lose structural identifiability when started from particular initial conditions (Denis-Vidal et al., 2001; Saccomani et al., 2003; Villaverde and Banga, 2017). Likewise, for the case of the inputs, structural identifiability analysis methods can determine whether a model is structurally identifiable provided that sufficiently exciting inputs are applied – but it is not straightforward to characterize the necessary inputs with the existing implementations of said methods. Our aim is to determine this qualitative information analytically, without performing a numerical optimal experimental design (OED) procedure. In this way we will be able to rule out insufficiently exciting inputs a priori, removing them from consideration in a subsequent OED. To the best of our knowledge no currently available tool informs about the shape of the input required for structural identifiability. In the present paper we fill this gap by building on an existing method, STRIKE-GOLDD (Villaverde et al., 2016). Here we show that, in its original form, STRIKE-GOLDD may classify an identifiable model as unidentifiable if a time-varying input is required for identification, which happens when a constant input does not excite the system dynamics sufficiently. We then modify the method by using extended Lie derivatives, which correctly analyses the effect of time-varying inputs. Importantly, in this way it is possible to inform about the type of input that is necessary in order to guarantee structural identifiability. We demonstrate this capability with a two-compartment model which is shown to require timevarying inputs in order for it to be structurally identifiable.
2.2 Structural identifiability as generalized observability Local structural identifiability can be considered as a generalization of observability. A model is observable if it is possible to determine its internal state x by observing its output y. Observability can be evaluated by calculating the rank of an observability matrix, OI , which represents a map between the model output y (and its derivatives y, ˙ y¨, . . .) on the one hand, and its state x on the other. For nonlinear models with constant input OI can be calculated using Lie derivatives, as follows:
L2f g(x) = ··· Lif g(x) =
∂g(x) f (x, u). ∂x
(4)
∂Lf g(x) f (x, u), ∂x ∂Li−1 f g(x) ∂x
(5)
f (x, u).
The observability rank condition (ORC) states that, if the system given by (1) with constant input u satisfies rank(O(x0 )) = n, with O defined by (3), then it is (locally) observable around x0 (Hermann and Krener, 1977).
(1)
The model’s structural identifiability can be evaluated in the same way as its observability. To this end, we consider the parameters pi as additional states with zero dynamics p˙ i = 0, i.e, we augment the state variable vector as x ˜ = (x, p). The augmented (or generalized) observabilityidentifiability matrix, OI (˜ x), is then defined as:
where f and g are analytic vector functions, p ∈ Rq is the parameter vector, u(t) ∈ Rr the input vector, x(t) ∈ Rn the state variable vector, and y(t) ∈ Rm the output vector. We drop the dependence on p for ease of notation. A parameter pi is structurally globally identifiable (s.g.i.) if it can be uniquely determined from the system output, that is, if for almost any p∗ ∈ Rq (i.e., for any p except those belonging to a set of measure zero) the following property holds for all t and all admissible inputs u (Ljung and Glad, 1994): y(t, pˆ) = y(t, p∗ ) ⇒ pˆi = p∗i
Subsequent derivatives can be recursively calculated:
Dynamic biological processes may be described by ordinary differential equations (ODEs) as follows: x(t) ˙ = f (x(t), u(t), p), y(t) = g(x(t), p), x0 = x(t0 , p)
∂ y(t) ∂x ∂ y(t) ˙ ∂x ∂ y¨(t) ∂x .. .
Lf g(x) =
2.1 Structural identifiability definitions
model M :
∂ g(x) ∂x ∂ (L g(x)) f ∂x ∂ 2 (3) O(x) = (Lf g(x)) = ∂x .. . ∂ (n−1) ∂ y (Ln−1 (t) g(x)) ∂x ∂x f where Lf g(x) is the Lie derivative of g with respect to f :
2. METHODOLOGY
17
x) = OI (˜
(2)
A parameter pi is structurally locally identifiable (s.l.i.) if for almost any p∗ there is a neighbourhood V (p∗ ) in which (2) holds. If (2) does not hold in any neighbourhood of p∗ , then pi is structurally unidentifiable (s.u.). A model M is s.g.i. if all its parameters are s.g.i.; it is s.u. if at least one of its parameters is s.u.; and it is s.l.i. if all its parameters are s.l.i. or s.g.i. and at least one of them is not s.g.i..
∂ g(˜ x) ∂x ˜ ∂ (Lf g(˜ x)) ∂x ˜ ∂ (L2 g(˜ x)) ∂x ˜ f .. . ∂ (Ln+q−1 g(˜ x)) ∂x ˜ f
(6)
The corresponding observability-identifiability condition (OIC) states that, if the system (1) with constant input u satisfies rank(OI (˜ x0 )) = n + q, with OI (˜ x0 ) given by (6), then it is locally observable and structurally locally identifiable in a neighbourhood N (˜ x0 ) of x ˜0 . 17
IFAC FOSBE 2018 18 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
2.3 Time-varying inputs and extended Lie derivatives
is sufficient for structural identifiability. If, however, OI is full rank for {u˙ = 0, u ¨ = 0}, but not for u˙ = 0, then a ramp input is necessary and sufficient, and a constant input is not. The effect of a specific input can be tested by entering the corresponding expression for u(t) in (6).
The input vector u(t) in (1) can in general be time-varying. However, the expressions for the Lie derivatives used in (4) and (5) implicitly assume that it is constant. For timevarying inputs, the definition of the Lie derivative must be modified in order to correspond to the output derivative. To this end we use an extended Lie derivative (Karlsson et al., 2012), which is defined by:
Lf g(˜ x) =
j=∞ ∂g(˜ ∂g(˜ x) x) (j+1) f (˜ x, u) + u . (j) ∂x ˜ ∂u j=0
3. EXAMPLE: TWO COMPARTMENT MODEL We illustrate the effect of using extended Lie derivatives with a biological model of a physiological system with two compartments (i.e. two states, of which one is measured) and one input. The model equations are given by:
(7)
j=∞ j=0
∂g(˜ x) (j+1) u =0 ∂u(j)
M1 :
Calculating the OI matrix of Eq. (6) for M1 with 5 extended Lie derivatives, as in (7–9), yields a 6 × 6 matrix with rank(OI ) = 6. Therefore the observabilityidentifiability rank condition (OIC) is satisfied. This means that the model is structurally identifiable (as long as the input is sufficiently exciting). If, however, the nonextended version of the Lie derivative is used (4–5), then rank(OI ) = 5, and the model is reported to be unidentifiable. This latter conclusion is correct only if the input is constant. The practical meaning of these results is that, in order to determine the values of the parameters (k1e , k12 , k21 , b) by measuring the model output (y = x1 ), it is necessary to perform an experiment with time-varying input (u˙ = 0). Moreover, we can characterize the type of time dependency that is needed for the input to enable structural identifiability. We can do this by replacing in OI higher order derivatives of the input with zero and recalculating rank(OI ). For the model M1 this procedure yields that rank(OI ) = 6 even if u ¨ = 0, as long as u˙ = 0; however, rank(OI ) reduces to 5 if u˙ = 0. Thus we know that a ramp input (u(t) = k1 · t + k2 ) suffices for structural identifiability of the parameters, but a constant one (u(t) = k) does not. Fig. 1 illustrates this fact.
(8)
and the extended Lie derivative defined above is identical to its non-extended counterpart. However, the summation term is not necessarily zero in higher order extended Lie derivatives, which are calculated as:
Lif g(˜ x) =
∂Li−1 x) f g(˜ ∂x ˜
f (˜ x, u) +
j=∞
∂Li−1 x) f g(˜
j=0
∂u(j)
u(j+1) (9)
This means that, if we calculate the matrix OI in (6) using extended Lie derivatives, we may encounter first order derivatives of the input u˙ in the third and subsequent rows (for one-dimensional outputs), second order derivatives u ¨ in the fourth and subsequent rows, and so on. Thus, when the structural identifiability analysis of a model requires two or more Lie derivatives, the result may indeed be affected by the presence of time-varying inputs. It should also be noted that, in practice, the calculation of j=∞ ∂Lif 1 g(˜x) (j+1) u can be truncated, since the ith Lie j=0 ∂u(j) derivative, Lif , depends on (u, u, ˙ . . . , u(i−2) ) but not on higher derivatives of the input. Hence, for i > 1 it suffices to calculate the extended Lie derivative as:
Lif g(˜ x) =
∂Li−1 x) f g(˜ ∂x ˜
f (˜ x, u) +
j=i−2 j=0
∂Li−1 x) f g(˜ ∂u(j)
x˙1 = −(k1e + k12 ) · x1 + k21 · x2 + b · u, x˙2 = k12 · x1 − k21 · x2 , (11) y = x1 where the unknown parameter vector is p = (k1e , k12 , k21 , b); the initial condition x2 (0) is also unknown. A diagram of the model is shown in the left hand side of Fig. 1.
In the second term of the sum within (7), u(j) and u(j+1) denote the j th and (j + 1)th derivatives of the input, respectively (note that we write u instead of u(t) to ease the notation). Clearly, since the output function g does not depend on the input directly (i.e. it is g(˜ x), not g(˜ x, u)), it holds that
4. CONCLUSION The method presented here analyses the structural identifiability of models with continuously time-varying inputs. It does so by including the derivatives of the input in the identifiability matrix through extended Lie derivatives. In this way, the approach takes into account the ability (or lack thereof) of a given time-varying input to excite a dynamic behaviour in the system that leads to the resolution of structural non-identifiabilities. This is important because, in the case of structural unidentifiability, parameter estimates are biologically meaningless, and model predictions may be wrong. The method can inform the design of new experiments: it delimits the type of external inputs that are required to correctly estimate the model parameters from the resulting dataset – more specifically, it establishes which derivatives of the input must be non-zero. It can also test the effect of a particular input by introducing its expression into the identifiability
u(j+1) (10)
2.4 Use of extended Lie derivatives for input design Using extended Lie derivatives (9, 10) we can characterize the type of time dependency that is needed for the input to enable structural identifiability. We can do this by setting to zero in OI (6) the derivatives of the input of order higher than a given one, and then recalculating rank(OI ). For example, if OI has full rank for u˙ = 0, a constant input 18
IFAC FOSBE 2018 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
19
y
u
b
k12 x1
x2 k21
k1e x1 = -(k1e + k12)·x1 + k21·x2 + b·u
x2 = k12·x1 - k21·x2 y = x1
Fig. 1. A two-compartment model. With a constant input it is not possible to distinguish between two different parameter vectors (note that the figure plots only two parameter vectors, but there is an infinite number of pairs of such vectors that are indistinguishable) and hence the parameters are structurally unidentifiable. In contrast, the parameters become identifiable with a time-varying input, e.g. with a ramp. This result does not depend on the number or location of the time-point measurements. matrix. We have demonstrated the methodology with a compartmental model that is structurally unidentifiable if a constant input is used, but becomes identifiable with a time-varying input such as a ramp function.
Ljung, L. and Glad, T. (1994). On global identifiability for arbitrary model parametrizations. Automatica, 30(2), 265–276. Meshkat, N., Kuo, C.E.z., and DiStefano III, J. (2014). On finding and using identifiable parameter combinations in nonlinear dynamic systems biology models and combos: A novel web implementation. PLoS One, 9(10). Pohjanpalo, H. (1978). System identifiability based on the power series expansion of the solution. Math. Biosci., 41(1), 21–33. Raue, A., Kreutz, C., Maiwald, T., Bachmann, J., et al. (2009). Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood. Bioinformatics, 25(15), 1923–1929. Raue, A., Steiert, B., Schelker, M., Kreutz, C., et al. (2015). Data2dynamics: a modeling environment tailored to parameter estimation in dynamical systems. Bioinformatics, 31(21), 3558–3560. Saccomani, M.P., Audoly, S., and D’Angi`o, L. (2003). Parameter identifiability of nonlinear systems: the role of initial conditions. Automatica, 39(4), 619–632. Saccomani, M.P. and Cobelli, C. (1992). Qualitative experiment design in physiological system identification. IEEE Control Syst. Mag., 12(6), 18–23. Stigter, J.D. and Molenaar, J. (2015). A fast algorithm to assess local structural identifiability. Automatica, 58, 118–124. Villaverde, A.F. and Banga, J.R. (2017). Structural properties of dynamic systems biology models: Identifiability, reachability, and initial conditions. Processes, 5(2), 29. Villaverde, A.F., Barreiro, A., and Papachristodoulou, A. (2016). Structural identifiability of dynamic systems biology models. PLoS Comput. Biol., 12(10), e1005153. ´ and Pronzato, L. (1990). Qualitative and Walter, E. quantitative experiment design for phenomenological models–a survey. Automatica, 26(2), 195–213. Xia, X. and Moog, C.H. (2003). Identifiability of nonlinear systems with application to HIV/AIDS models. IEEE Trans. Autom. Control, 48(2), 330–336. Yates, J.W., Evans, N.D., and Chappell, M.J. (2009). Structural identifiability analysis via symmetries of differential equations. Automatica, 45(11), 2585–2591.
REFERENCES Anguelova, M., Karlsson, J., and Jirstrand, M. (2012). Minimal output sets for identifiability. Math. Biosci., 239(1), 139–153. August, E. (2009). Parameter identifiability and optimal experimental design. In Int. Conf. on Computational Science and Engineering, CSE’09., volume 1, 277–284. Bellu, G., Saccomani, M.P., Audoly, S., and D´ Angio, L. (2007). DAISY: a new software tool to test global identifiability of biological and physiological systems. Comput. Methods Programs Biomed., 88(1), 52–61. Chi¸s, O., Banga, J.R., and Balsa-Canto, E. (2011). GenSSI: a software toolbox for structural identifiability analysis of biological models. Bioinformatics, 27(18), 2610–2611. Cobelli, C. and DiStefano, J. (1980). Parameter and structural identifiability concepts and ambiguities: a critical review and analysis. Am. J. Physiol. Regul. Integr. Comp. Physiol., 239(1), R7–R24. Denis-Vidal, L. and Joly-Blanchard, G. (2000). An easy to check criterion for (un)identifiability of uncontrolled systems and its applications. IEEE Trans. Autom. Control, 45, 768–771. Denis-Vidal, L., Joly-Blanchard, G., and Noiret, C. (2001). Some effective approaches to check the identifiability of uncontrolled nonlinear systems. Math. Comput. Simul, 57(1), 35–44. Evans, N.D., Chapman, M.J., Chappell, M.J., and Godfrey, K.R. (2002). Identifiability of uncontrolled nonlinear rational systems. Automatica, 38(10), 1799–1805. Hermann, R. and Krener, A.J. (1977). Nonlinear controllability and observability. IEEE Trans. Autom. Control, 22(5), 728–740. Karlsson, J., Anguelova, M., and Jirstrand, M. (2012). An efficient method for structural identiability analysis of large dynamic systems. In 16th IFAC Symposium on System Identification, volume 16, 941–946. 19
ScienceDirect
IFAC PapersOnLine 51-19 (2018) 16–19
Sufficiently Exciting Inputs for Structurally Sufficiently Exciting Inputs for Structurally Sufficiently Exciting Inputs for Structurally Identifiable Systems Biology Models Sufficiently Exciting Inputs for Structurally Identifiable Systems Biology Identifiable Systems Biology Models Models Identifiable Systems Biology Models ∗ Alejandro F. Villaverde Neil D. Evans ∗∗
∗ ∗∗ Alejandro F. Villaverde ∗∗ ∗ Neil D. Evans∗∗∗ Michael J. Julio Alejandro F.Chappell Villaverde NeilR. D.Banga Evans∗∗∗ ∗∗ ∗ Michael J. Chappell Julio R. Banga Alejandro Villaverde NeilR. D.Banga Evans∗ ∗∗ Julio Michael J.F.Chappell ∗∗ ∗ Michael J. Chappell Julio R. Banga ∗ ∗ Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain ∗ Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain Bioprocess Engineering Group, IIM-CSIC, Vigo 36208, Spain (e-mail: ∗ (e-mail: [email protected]). [email protected]). ∗∗ Bioprocess Engineering Group, IIM-CSIC, 36208,CV4 Spain University of Warwick,Vigo Coventry 7AL, (e-mail: [email protected]). ∗∗ School of Engineering, University of Warwick, Coventry CV4 7AL, (e-mail: [email protected]). ∗∗ School of Engineering, UK School of Engineering, University of Warwick, Coventry CV4 ∗∗ UK of Warwick, Coventry CV4 7AL, School of Engineering, University 7AL, UK UK Abstract: its value can theoretically be estimated Abstract: A A parameter parameter is is structurally structurally identifiable identifiable if if its value can theoretically be estimated by observing the model output. Structural identifiability is a desirable property be in estimated biological Abstract: A parameter is structurally identifiable if its value can theoretically by observingA the model output. Structural identifiability is aa desirable property in biological Abstract: parameter is structurally identifiable if its value can theoretically be by observingif the model output. Structural identifiability is its desirable property in estimated biological modelling: a parameter is structurally unidentifiable, estimated numerical value modelling: if a parameter is structurally unidentifiable, its estimated numerical value is is by observingifand the model predictions output. Structural identifiability is its a desirable property in biological modelling: a parameter is structurally unidentifiable, estimated numerical value is meaningless, model of unmeasured state variables can be wrong, compromising meaningless, and model predictions of unmeasured state variables can be wrong, compromising modelling: if a parameter is structurally unidentifiable, its estimated numerical value is the ability the model to insight. Structural on the meaningless, of unmeasured state variablesidentifiability can be wrong,depends compromising the ability of of and the model model predictions to provide provide biological biological insight. Structural identifiability depends on the meaningless, and model predictions of unmeasured state variables can be wrong, compromising system dynamics, observation function (model output), initial and the ability of the model to provide biological Structural identifiability depends inputs. on the system dynamics, observation function (modelinsight. output), initial conditions, conditions, and external external inputs. the ability of we thefocus model to provide biological insight. Structural identifiability depends on the system dynamics, observation function (model output), initial conditions, andanalysis external inputs. In this paper on the last factor. Methods for structural identifiability typically In this paper we focus on the last factor.(model Methods for structural identifiability analysis typically system dynamics, observation function output), initial conditions, and external inputs. classify a model as identifiable provided that it is fed with sufficiently exciting inputs. For In this paper we focus on the last factor. Methods for structural identifiability analysis typically classify a model as identifiable provided that it for is structural fed with sufficiently For exciting inputs. In this paper we focus onmay the last factor. Methods identifiability analysis typically example, a given model require a time-varying input to be structurally identifiable, while classify a model as identifiable provided that it is fed with sufficiently exciting inputs. For example,a amodel given model may require a time-varying input to besufficiently structurally identifiable, while classify as identifiable provided that it is fed with exciting inputs. For for another model a constant non-zero input may be enough. Here we present a method that example, a given model may require a time-varying input to be structurally identifiable, while for another modelmodel a constant non-zero input may beinput enough. Here we present aa method that example, a given may require a time-varying to be structurally identifiable, while for another model a constant non-zero input may be enough. Here we present method that determines how sufficiently exciting an should be. The builds on STRIKEdetermines how sufficiently exciting an input input The approach approach onathe the STRIKEfor another model a constant non-zero input should may bebe. enough. Here we builds present method that GOLDD toolbox, which considers structural identifiability as generalized observability. The determines how sufficiently exciting an input should be. The approach builds on the STRIKEGOLDD toolbox, which considers structural identifiability as generalized observability. The determines how sufficiently exciting an input should be. The approach builds on the STRIKEapproach incorporates extended Lie derivatives, which correctly assess structural identifiability GOLDD toolbox, which considers structural identifiability as generalized observability. The approach incorporates extended Lie structural derivatives,identifiability which correctly assess structural identifiability GOLDD toolbox, which considers as generalized observability. The in the case of time-varying inputs. The procedure can also be used to determine the type approach incorporates extended Lie derivatives, which correctly assess structural identifiability in the case of time-varying inputs. The procedure cancorrectly also be used tostructural determineidentifiability the type of of approach incorporates extended Lie derivatives, which assess in the case of time-varying inputs. The procedure can also be used to determine the type of input profile that is to the identifiable. This capability is helpful when input profile is required required inputs. to make makeThe the parameters parameters identifiable. that This capability is helpful when in the case of time-varying procedure can also be used to determine the type of designing new experiments for the purpose of parameter estimation. input profile that is required to make the parameters identifiable. This capability is helpful when designing new experiments for the purpose of parameter estimation. input profile that is required to make the parameters identifiable. This capability is helpful when designing new experiments for the purpose of parameter estimation. © 2018, IFAC of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. designing new(International experimentsFederation for the purpose of parameter estimation. Keywords: properties, Observability, Keywords: Identifiability, Identifiability, Structural Structural Identifiability, Identifiability, Structural Structural properties, Observability, Keywords: Identifiability, Structural Identifiability, Structural properties, Observability, Parameter identification, estimation, Parameter identification, Parameter Parameter estimation, Inputs. Inputs. Keywords: Identifiability, Structural Identifiability, Structural properties, Observability, Parameter identification, Parameter estimation, Inputs. Parameter identification, Parameter estimation, Inputs. 1. INTRODUCTION Molenaar, 2015). There are a number of software tools 1. tools Molenaar, There are aa number software 1. INTRODUCTION INTRODUCTION Molenaar, 2015). 2015). There are number of of methodologies, software tools implementing some of the aforementioned implementing some of the aforementioned methodologies, 1. INTRODUCTION Molenaar, 2015). There are a number of software tools A parameter is structurally identifiable if it is theoreti- implementing such as DAISY (Bellu et al., 2007), GenSSI (Chi¸ ss et al., some of the aforementioned methodologies, A parameter is structurally identifiable if it is theoretisuch as DAISY (Bellu et al., 2007), GenSSI (Chi¸ et al., implementing some of the aforementioned methodologies, A parameter is structurally identifiable if it is theoretically possible to determine its true value from noiseless such 2011), EAR (Karlsson 2012), COMBOS (Meshkat as DAISY (Bellu et al., 2007), GenSSI (Chi¸ s et al., cally possible to determine its true value from noiseless 2011), EAR (Karlsson et al., 2012), COMBOS (Meshkat A parameter is structurally identifiable if it is theoretias2014), DAISY (Bellu et al., 2007), GenSSI (Chi¸ s 2016), et al., cally possible and to determine true value noiseless such data (Cobelli DiStefano,its1980). Since from identifiability 2011), EAR (Karlsson et al., 2012), COMBOS (Meshkat et al., STRIKE-GOLDD (Villaverde et al., data and DiStefano, 1980). Since identifiability et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), cally (Cobelli possible to other determine true noiseless EAR (Karlsson et al.,et2012), COMBOS (Meshkat data (Cobelli and DiStefano, 1980). Since from identifiability depends, among things,itson the value system’s inputs and 2011), and Data2Dynamics (Raue al., 2015). These methods et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), depends, among other things, on the system’s inputs and and Data2Dynamics (Raue et al., 2015). These methods data (Cobelli and DiStefano, 1980). Since identifiability et al., 2014), STRIKE-GOLDD (Villaverde et al., 2016), depends, among other things, on the system’s inputs and outputs, it can vary as a result of the experiment design. and assess whether a given input-output configuration is in Data2Dynamics (Raue et al., 2015). These methods outputs, it can as a the experiment design. whether a given input-output configuration is in depends, other things, the and assess and Data2Dynamics (Raue et al., 2015). These methods outputs, can vary vary as a result resultonof of thesystem’s experiment The termitamong qualitative experiment design refers toinputs thedesign. selecprinciple sufficient to obtain a structurally identifiable assess whether a given input-output configuration is in The term qualitative experiment design refers to the selecprinciple sufficient to obtain a structurally identifiable outputs, it can vary as a result of the experiment design. assess whether a given input-output configuration is in The qualitative experiment design to the selection term of input and output ports in orderrefers to maximize the principle sufficient to Cobelli, obtain a1992). structurally identifiable model (Saccomani and Algorithms that tion of input and output ports in order to maximize the model (Saccomani and Cobelli, Algorithms that auauThe term experiment design to the selecsufficient tothe obtain a1992). structurally identifiable tion of input and output ports in orderrefers to maximize the principle number ofqualitative structurally identifiable parameters (Walter and tomate the search for subsets of outputs that guarantee model (Saccomani and Cobelli, 1992). Algorithms that aunumber of structurally identifiable parameters (Walter and tomate the search for the subsets of outputs that guarantee tion of input and output ports in order to maximize the model (Saccomani and Cobelli, 1992). Algorithms that aunumber of structurally identifiable parameters (Walter and Pronzato, 1990). structural identifiability are also available (August, 2009; tomate the search for the subsets of outputs that guarantee Pronzato, identifiability are also ofavailable (August, 2009; number of 1990). structurally identifiable parameters (Walter and structural tomate the search for the subsets outputs that guarantee Pronzato, 1990). Anguelova et identifiability Some structural identifiability methods, such as the simi- structural Anguelova et al., al., 2012). 2012). are also available (August, 2009; Pronzato, 1990). identifiability structural Some structural methods, such as the simiAnguelova identifiability et al., 2012). are also available (August, 2009; Some structural identifiability methods, such as the similarity transformation approach (Evans et al., 2002; Yates Anguelova Determining etthe al.,necessary 2012). set larity transformation approach (Evans et al., 2002; the necessary set of of inputs inputs and and outputs outputs is is only only Some identifiability the Yates simi- Determining larity transformation approach (Evans etsuch al., as 2002; Yates et al.,structural 2009) and direct testmethods, (Denis-Vidal and Jolypart of experiment design. Additionally, it is necessary Determining the necessary set of inputs and outputs is only et al., 2009) and direct test (Denis-Vidal and Jolypart of experiment design. Additionally, it is necessary larity transformation approach (Evans et al., 2002; Yates Determining the necessary setAdditionally, ofnumber inputs and it outputs is only et al., 2009) andaredirect test (Denis-Vidal Joly- part Blanchard, 2000) applicable to autonomousand systems, of experiment design. is necessary to aspects as of Blanchard, 2000) aredirect applicable to systems, to define define aspects such such as the theAdditionally, number and andittiming timing of the the et systems al., 2009) and test approaches, (Denis-Vidal and Joly- part of experiment design. is necessary Blanchard, 2000) applicable to autonomous autonomous systems, i.e. with noare input. Other such as those samples, the system’s initial state, and characterization of to define aspects such as the number and timing of the i.e. systems with no input. Other approaches, such as those samples, the system’s initial state, and characterization of Blanchard, 2000) are applicable to autonomous systems, to define aspects such as the number and timing of the i.e. systems with no input. Other approaches, such as those based on power series (Pohjanpalo, 1978), differential al- samples, the inputs. This task is known as quantitative experiment the system’s initial state, and characterization of based on series (Pohjanpalo, 1978), differential althe inputs. This task isinitial known as quantitative experiment i.e. systems with noGlad, input. Other approaches, such as those samples, the system’s state, and characterization of based on power power (Pohjanpalo, 1978), differential al- design, gebra (Ljung andseries 1994; Bellu et al., 2007), implicit and aims at maximizing practical identifiability. the inputs. This task is known as quantitative experiment gebra andseries Glad,(Pohjanpalo, 1994; Bellu et al., implicit andThis aimstask at maximizing practical identifiability. based (Ljung on power differential al- design, the the inputs. is we known as quantitative experiment gebra (Ljung 1994; Bellu et al., 2007), 2007), implicit functions (Xiaand andGlad, Moog, 2003), or 1978), differential geometry design, and aims at maximizing practical identifiability. In present paper address a different, although functions (Xia and Moog, 2003), or differential geometry In the present we addresspractical a different, although gebra (Ljung 1994; Bellu al.,to2007), implicit and aimspaper at maximizing identifiability. functions (Xia andGlad, Moog, 2003), or et differential geometry (Villaverde etand al., 2016), can be applied systems with design, related, problem. We seek to determine, a priori, what In the present paper we address a different, although (Villaverde et al., 2016), can be applied to systems with related, problem. We seek to determine, a priori, what functions (Xia and Moog, 2003), or differential geometry In the present paper we address a different, although (Villaverde et al., 2016), can be applied to systems with external inputs. Additionally, it is possible (for a given set related, mathematical form the inputs must have in order to guarproblem. We seek to determine, a priori, what external inputs. Additionally, it is possible (for a given set mathematical form the inputs must have in order to guar(Villaverde et al., 2016), can be applied to systems with problem. We seek to must determine, aorder priori, what external inputs.with Additionally, it is possible a given set related, of experiments defined inputs) to assess(for identifiability mathematical form the inputs have in to guarantee structural identifiability. For example, for a given of experiments defined to assess(for identifiability structural identifiability. For example, for a given external inputs.with Additionally, it is possible a given set antee mathematical form the inputsconstant must have in order to guarof experiments with defined inputs) inputs) to identifiability with numerical approaches based onassess profile likelihoods antee structural identifiability. For example, for a given model structure a non-zero input may not be with numerical approaches based on profile likelihoods model structural structure identifiability. a non-zero constant input may be of experiments with defined inputs) assess identifiability For example, for anot given with numerical approaches based to onmatrix profile likelihoods (Raue et al., 2009) or the sensitivity (Stigter and antee sufficiently exciting and a ramp input may be necessary. To model structure a non-zero constant input may not be (Raue et al., 2009) or the sensitivity and sufficiently exciting and a ramp input may be necessary. To with numerical approaches based onmatrix profile(Stigter likelihoods model structure a non-zero constant input may not be (Raue et al., 2009) or the sensitivity matrix (Stigter and aa certain extent, this situation resembles the relationship This project has received funding from the European Union’s Horisufficiently exciting and a ramp input may be necessary. To certain extent, this situation resembles thenecessary. relationship (Raue et al.,has 2009) or the sensitivity matrix (Stigter and sufficiently exciting and a ramp input may be To This project received funding from the European Union’s Horia certain extent, this situation resembles the relationship between structural identifiability and conditions: zon 2020 research innovation programme under grant agreement This project hasand received funding from the European Union’s Horibetween identifiability and initial initial conditions: zon 2020 research innovation programme under grant agreement a certain structural extent, this situation resembles theidentifiability relationship This project hasand received funding from the European Union’s Horigenerally, methods that analyse structural between structural identifiability and initial conditions: No 686282 (“CANPATHPRO”). zon 2020 research and innovation programme under grant agreement generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). between structural identifiability and initial conditions: zon 2020 research and innovation programme under grant agreement generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). generally, methods that analyse structural identifiability No 686282 (“CANPATHPRO”). 2405-8963 © 2018, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved.
Copyright © 2018 IFAC 16 Peer review©under responsibility of International Federation of Automatic Copyright 2018 IFAC 16 Control. Copyright © 2018 IFAC 16 10.1016/j.ifacol.2018.09.015 Copyright © 2018 IFAC 16
IFAC FOSBE 2018 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
yield results that are valid for almost all initial conditions; however, a model classified as structurally identifiable may lose structural identifiability when started from particular initial conditions (Denis-Vidal et al., 2001; Saccomani et al., 2003; Villaverde and Banga, 2017). Likewise, for the case of the inputs, structural identifiability analysis methods can determine whether a model is structurally identifiable provided that sufficiently exciting inputs are applied – but it is not straightforward to characterize the necessary inputs with the existing implementations of said methods. Our aim is to determine this qualitative information analytically, without performing a numerical optimal experimental design (OED) procedure. In this way we will be able to rule out insufficiently exciting inputs a priori, removing them from consideration in a subsequent OED. To the best of our knowledge no currently available tool informs about the shape of the input required for structural identifiability. In the present paper we fill this gap by building on an existing method, STRIKE-GOLDD (Villaverde et al., 2016). Here we show that, in its original form, STRIKE-GOLDD may classify an identifiable model as unidentifiable if a time-varying input is required for identification, which happens when a constant input does not excite the system dynamics sufficiently. We then modify the method by using extended Lie derivatives, which correctly analyses the effect of time-varying inputs. Importantly, in this way it is possible to inform about the type of input that is necessary in order to guarantee structural identifiability. We demonstrate this capability with a two-compartment model which is shown to require timevarying inputs in order for it to be structurally identifiable.
2.2 Structural identifiability as generalized observability Local structural identifiability can be considered as a generalization of observability. A model is observable if it is possible to determine its internal state x by observing its output y. Observability can be evaluated by calculating the rank of an observability matrix, OI , which represents a map between the model output y (and its derivatives y, ˙ y¨, . . .) on the one hand, and its state x on the other. For nonlinear models with constant input OI can be calculated using Lie derivatives, as follows:
L2f g(x) = ··· Lif g(x) =
∂g(x) f (x, u). ∂x
(4)
∂Lf g(x) f (x, u), ∂x ∂Li−1 f g(x) ∂x
(5)
f (x, u).
The observability rank condition (ORC) states that, if the system given by (1) with constant input u satisfies rank(O(x0 )) = n, with O defined by (3), then it is (locally) observable around x0 (Hermann and Krener, 1977).
(1)
The model’s structural identifiability can be evaluated in the same way as its observability. To this end, we consider the parameters pi as additional states with zero dynamics p˙ i = 0, i.e, we augment the state variable vector as x ˜ = (x, p). The augmented (or generalized) observabilityidentifiability matrix, OI (˜ x), is then defined as:
where f and g are analytic vector functions, p ∈ Rq is the parameter vector, u(t) ∈ Rr the input vector, x(t) ∈ Rn the state variable vector, and y(t) ∈ Rm the output vector. We drop the dependence on p for ease of notation. A parameter pi is structurally globally identifiable (s.g.i.) if it can be uniquely determined from the system output, that is, if for almost any p∗ ∈ Rq (i.e., for any p except those belonging to a set of measure zero) the following property holds for all t and all admissible inputs u (Ljung and Glad, 1994): y(t, pˆ) = y(t, p∗ ) ⇒ pˆi = p∗i
Subsequent derivatives can be recursively calculated:
Dynamic biological processes may be described by ordinary differential equations (ODEs) as follows: x(t) ˙ = f (x(t), u(t), p), y(t) = g(x(t), p), x0 = x(t0 , p)
∂ y(t) ∂x ∂ y(t) ˙ ∂x ∂ y¨(t) ∂x .. .
Lf g(x) =
2.1 Structural identifiability definitions
model M :
∂ g(x) ∂x ∂ (L g(x)) f ∂x ∂ 2 (3) O(x) = (Lf g(x)) = ∂x .. . ∂ (n−1) ∂ y (Ln−1 (t) g(x)) ∂x ∂x f where Lf g(x) is the Lie derivative of g with respect to f :
2. METHODOLOGY
17
x) = OI (˜
(2)
A parameter pi is structurally locally identifiable (s.l.i.) if for almost any p∗ there is a neighbourhood V (p∗ ) in which (2) holds. If (2) does not hold in any neighbourhood of p∗ , then pi is structurally unidentifiable (s.u.). A model M is s.g.i. if all its parameters are s.g.i.; it is s.u. if at least one of its parameters is s.u.; and it is s.l.i. if all its parameters are s.l.i. or s.g.i. and at least one of them is not s.g.i..
∂ g(˜ x) ∂x ˜ ∂ (Lf g(˜ x)) ∂x ˜ ∂ (L2 g(˜ x)) ∂x ˜ f .. . ∂ (Ln+q−1 g(˜ x)) ∂x ˜ f
(6)
The corresponding observability-identifiability condition (OIC) states that, if the system (1) with constant input u satisfies rank(OI (˜ x0 )) = n + q, with OI (˜ x0 ) given by (6), then it is locally observable and structurally locally identifiable in a neighbourhood N (˜ x0 ) of x ˜0 . 17
IFAC FOSBE 2018 18 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
2.3 Time-varying inputs and extended Lie derivatives
is sufficient for structural identifiability. If, however, OI is full rank for {u˙ = 0, u ¨ = 0}, but not for u˙ = 0, then a ramp input is necessary and sufficient, and a constant input is not. The effect of a specific input can be tested by entering the corresponding expression for u(t) in (6).
The input vector u(t) in (1) can in general be time-varying. However, the expressions for the Lie derivatives used in (4) and (5) implicitly assume that it is constant. For timevarying inputs, the definition of the Lie derivative must be modified in order to correspond to the output derivative. To this end we use an extended Lie derivative (Karlsson et al., 2012), which is defined by:
Lf g(˜ x) =
j=∞ ∂g(˜ ∂g(˜ x) x) (j+1) f (˜ x, u) + u . (j) ∂x ˜ ∂u j=0
3. EXAMPLE: TWO COMPARTMENT MODEL We illustrate the effect of using extended Lie derivatives with a biological model of a physiological system with two compartments (i.e. two states, of which one is measured) and one input. The model equations are given by:
(7)
j=∞ j=0
∂g(˜ x) (j+1) u =0 ∂u(j)
M1 :
Calculating the OI matrix of Eq. (6) for M1 with 5 extended Lie derivatives, as in (7–9), yields a 6 × 6 matrix with rank(OI ) = 6. Therefore the observabilityidentifiability rank condition (OIC) is satisfied. This means that the model is structurally identifiable (as long as the input is sufficiently exciting). If, however, the nonextended version of the Lie derivative is used (4–5), then rank(OI ) = 5, and the model is reported to be unidentifiable. This latter conclusion is correct only if the input is constant. The practical meaning of these results is that, in order to determine the values of the parameters (k1e , k12 , k21 , b) by measuring the model output (y = x1 ), it is necessary to perform an experiment with time-varying input (u˙ = 0). Moreover, we can characterize the type of time dependency that is needed for the input to enable structural identifiability. We can do this by replacing in OI higher order derivatives of the input with zero and recalculating rank(OI ). For the model M1 this procedure yields that rank(OI ) = 6 even if u ¨ = 0, as long as u˙ = 0; however, rank(OI ) reduces to 5 if u˙ = 0. Thus we know that a ramp input (u(t) = k1 · t + k2 ) suffices for structural identifiability of the parameters, but a constant one (u(t) = k) does not. Fig. 1 illustrates this fact.
(8)
and the extended Lie derivative defined above is identical to its non-extended counterpart. However, the summation term is not necessarily zero in higher order extended Lie derivatives, which are calculated as:
Lif g(˜ x) =
∂Li−1 x) f g(˜ ∂x ˜
f (˜ x, u) +
j=∞
∂Li−1 x) f g(˜
j=0
∂u(j)
u(j+1) (9)
This means that, if we calculate the matrix OI in (6) using extended Lie derivatives, we may encounter first order derivatives of the input u˙ in the third and subsequent rows (for one-dimensional outputs), second order derivatives u ¨ in the fourth and subsequent rows, and so on. Thus, when the structural identifiability analysis of a model requires two or more Lie derivatives, the result may indeed be affected by the presence of time-varying inputs. It should also be noted that, in practice, the calculation of j=∞ ∂Lif 1 g(˜x) (j+1) u can be truncated, since the ith Lie j=0 ∂u(j) derivative, Lif , depends on (u, u, ˙ . . . , u(i−2) ) but not on higher derivatives of the input. Hence, for i > 1 it suffices to calculate the extended Lie derivative as:
Lif g(˜ x) =
∂Li−1 x) f g(˜ ∂x ˜
f (˜ x, u) +
j=i−2 j=0
∂Li−1 x) f g(˜ ∂u(j)
x˙1 = −(k1e + k12 ) · x1 + k21 · x2 + b · u, x˙2 = k12 · x1 − k21 · x2 , (11) y = x1 where the unknown parameter vector is p = (k1e , k12 , k21 , b); the initial condition x2 (0) is also unknown. A diagram of the model is shown in the left hand side of Fig. 1.
In the second term of the sum within (7), u(j) and u(j+1) denote the j th and (j + 1)th derivatives of the input, respectively (note that we write u instead of u(t) to ease the notation). Clearly, since the output function g does not depend on the input directly (i.e. it is g(˜ x), not g(˜ x, u)), it holds that
4. CONCLUSION The method presented here analyses the structural identifiability of models with continuously time-varying inputs. It does so by including the derivatives of the input in the identifiability matrix through extended Lie derivatives. In this way, the approach takes into account the ability (or lack thereof) of a given time-varying input to excite a dynamic behaviour in the system that leads to the resolution of structural non-identifiabilities. This is important because, in the case of structural unidentifiability, parameter estimates are biologically meaningless, and model predictions may be wrong. The method can inform the design of new experiments: it delimits the type of external inputs that are required to correctly estimate the model parameters from the resulting dataset – more specifically, it establishes which derivatives of the input must be non-zero. It can also test the effect of a particular input by introducing its expression into the identifiability
u(j+1) (10)
2.4 Use of extended Lie derivatives for input design Using extended Lie derivatives (9, 10) we can characterize the type of time dependency that is needed for the input to enable structural identifiability. We can do this by setting to zero in OI (6) the derivatives of the input of order higher than a given one, and then recalculating rank(OI ). For example, if OI has full rank for u˙ = 0, a constant input 18
IFAC FOSBE 2018 Chicago, Illinois, USA, August 5-8, 2018 Alejandro F. Villaverde et al. / IFAC PapersOnLine 51-19 (2018) 16–19
19
y
u
b
k12 x1
x2 k21
k1e x1 = -(k1e + k12)·x1 + k21·x2 + b·u
x2 = k12·x1 - k21·x2 y = x1
Fig. 1. A two-compartment model. With a constant input it is not possible to distinguish between two different parameter vectors (note that the figure plots only two parameter vectors, but there is an infinite number of pairs of such vectors that are indistinguishable) and hence the parameters are structurally unidentifiable. In contrast, the parameters become identifiable with a time-varying input, e.g. with a ramp. This result does not depend on the number or location of the time-point measurements. matrix. We have demonstrated the methodology with a compartmental model that is structurally unidentifiable if a constant input is used, but becomes identifiable with a time-varying input such as a ramp function.
Ljung, L. and Glad, T. (1994). On global identifiability for arbitrary model parametrizations. Automatica, 30(2), 265–276. Meshkat, N., Kuo, C.E.z., and DiStefano III, J. (2014). On finding and using identifiable parameter combinations in nonlinear dynamic systems biology models and combos: A novel web implementation. PLoS One, 9(10). Pohjanpalo, H. (1978). System identifiability based on the power series expansion of the solution. Math. Biosci., 41(1), 21–33. Raue, A., Kreutz, C., Maiwald, T., Bachmann, J., et al. (2009). Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood. Bioinformatics, 25(15), 1923–1929. Raue, A., Steiert, B., Schelker, M., Kreutz, C., et al. (2015). Data2dynamics: a modeling environment tailored to parameter estimation in dynamical systems. Bioinformatics, 31(21), 3558–3560. Saccomani, M.P., Audoly, S., and D’Angi`o, L. (2003). Parameter identifiability of nonlinear systems: the role of initial conditions. Automatica, 39(4), 619–632. Saccomani, M.P. and Cobelli, C. (1992). Qualitative experiment design in physiological system identification. IEEE Control Syst. Mag., 12(6), 18–23. Stigter, J.D. and Molenaar, J. (2015). A fast algorithm to assess local structural identifiability. Automatica, 58, 118–124. Villaverde, A.F. and Banga, J.R. (2017). Structural properties of dynamic systems biology models: Identifiability, reachability, and initial conditions. Processes, 5(2), 29. Villaverde, A.F., Barreiro, A., and Papachristodoulou, A. (2016). Structural identifiability of dynamic systems biology models. PLoS Comput. Biol., 12(10), e1005153. ´ and Pronzato, L. (1990). Qualitative and Walter, E. quantitative experiment design for phenomenological models–a survey. Automatica, 26(2), 195–213. Xia, X. and Moog, C.H. (2003). Identifiability of nonlinear systems with application to HIV/AIDS models. IEEE Trans. Autom. Control, 48(2), 330–336. Yates, J.W., Evans, N.D., and Chappell, M.J. (2009). Structural identifiability analysis via symmetries of differential equations. Automatica, 45(11), 2585–2591.
REFERENCES Anguelova, M., Karlsson, J., and Jirstrand, M. (2012). Minimal output sets for identifiability. Math. Biosci., 239(1), 139–153. August, E. (2009). Parameter identifiability and optimal experimental design. In Int. Conf. on Computational Science and Engineering, CSE’09., volume 1, 277–284. Bellu, G., Saccomani, M.P., Audoly, S., and D´ Angio, L. (2007). DAISY: a new software tool to test global identifiability of biological and physiological systems. Comput. Methods Programs Biomed., 88(1), 52–61. Chi¸s, O., Banga, J.R., and Balsa-Canto, E. (2011). GenSSI: a software toolbox for structural identifiability analysis of biological models. Bioinformatics, 27(18), 2610–2611. Cobelli, C. and DiStefano, J. (1980). Parameter and structural identifiability concepts and ambiguities: a critical review and analysis. Am. J. Physiol. Regul. Integr. Comp. Physiol., 239(1), R7–R24. Denis-Vidal, L. and Joly-Blanchard, G. (2000). An easy to check criterion for (un)identifiability of uncontrolled systems and its applications. IEEE Trans. Autom. Control, 45, 768–771. Denis-Vidal, L., Joly-Blanchard, G., and Noiret, C. (2001). Some effective approaches to check the identifiability of uncontrolled nonlinear systems. Math. Comput. Simul, 57(1), 35–44. Evans, N.D., Chapman, M.J., Chappell, M.J., and Godfrey, K.R. (2002). Identifiability of uncontrolled nonlinear rational systems. Automatica, 38(10), 1799–1805. Hermann, R. and Krener, A.J. (1977). Nonlinear controllability and observability. IEEE Trans. Autom. Control, 22(5), 728–740. Karlsson, J., Anguelova, M., and Jirstrand, M. (2012). An efficient method for structural identiability analysis of large dynamic systems. In 16th IFAC Symposium on System Identification, volume 16, 941–946. 19