Copyright © IF AC Modelling and Control in Biomedical Systems, Warwick, UK, 1997
REPARAMETERISATION OF UNIDENTIFIABLE SYSTEMS USING THE TAYLOR SERIES APPROACH Roger N. Gunn· Michael J. Chappell·· Vincent J. Cunningham·
• Cyclotron Unit, MRC Clinical Sciences Centre, Royal Postgraduate Medical School, Hammersmith Hospital, Du Cane Road, London, W12 ONN, England . •• Department of Engineering, University of Warwick, Coventry, CV4 7AL, England.
Abstract: This paper presents theory for the e."Cistence of a minimal locally identifiable reparameterisation of an unidentifiable system and is based on the Taylor series. approach of Pohjanpalo. Furthermore, a methodology is given which facilitates the construction of a reparameterisation. The reparameterisation reduces a system to its minimal form (in the sense of the number of parameters in the reparameterisation) and ensures that all the parameters are (at least) locally identifiable. The method for determining the parameter vector dimension of the minimal systems involves extending the Jacobian rank test of Pohjanpalo. Keywords: Identifiability, Models, Parameters, Parametrization, Methodology, Algebraic Approaches
yielded by the transfer function . The rank of the system's Jacobian matri..x can be utilised in assessing the identifiability of the system parameters (Rothenburg, 1971) . Parameter bounds can be derived for linear systems by using constraints implicit in a compartmental structure (Cobelli and Toffolo, 1987). The e.'Chaustive modelling approach has been used to obtain the whole set of output indistinguishable linear compartmental models (WaIter and Lecourtier, 1981).
1. INTRODUCTION
The determination of the structural identifiability of system parameters is an important procedure since it is a necessary prerequisite for experimental design. Over recent years a large amount of work has centred on the development of analytical methods and approaches for the determination of the identifiability of system parameters both for linear and nonlinear models (see for example (Godfrey and DiStefano, 1987) (Pohjanpalo, 1978)) . An important question is ' what can be done with unidentifiable parameters and systems ?' . Aside from obvious model simplification or redesign of the experiment considered, little has been conjectured theoretically. A method for analysing unidentifiable linear models which involves classifying the parameters as identifiable or having some degree of freedom has been presented (Milanese and Sorrentino, 1978). This is achieved by looking at the coefficients
In this paper, a method is presented which offers criteria for the e."Cistence of a locally identifiable reparameterised version of an unidentifiable system. A method for determining the reparameterisation using locally identifiable parameter combinations and a state space transformation is introduced . This reparameterised system is essentially a reduction of the original system to its minimal locally identifiable form (with respect to the parameterisation). The Taylor series approach (Pohjanpalo, 1978) and Pohjanpalo's Jacobian
247
rank test (Pohjanpalo, 1982) are employed which means that the approach is applicable to both linear and nonlinear systems and also systems with specified inputs. The technique relies on the calculation of the rank of the system 's (possibly infinite) J acobian matrix, and if this is intractable the procedure cannot be applied. The whole procedure is envisaged to be fully applied only when the new (at least locally identifiable) parameters are relevant to the underlying physical system, otherwise redesign of the experiment is a requisite. The examples presented in this paper will only consider systems where the appropriate Jacobian matrix is rank deficient by 1, although the results presented here apply for higher dimensions of rank deficiency. Illustrative examples of both linear and nonlinear biological systems are given. 2. BACKGROUND THEORY
parameterised system assuming that y and its derivatives, evaluated at some t = to, are available and perfect (noiseless) . These are essentially measurable constants (and in particular are independent of B). These coefficients are denoted by
[ y~:Ol] ,[Y~:ll] ,..., [y\nl] ,. .. y~)
y~)
(2)
y~)
and the Taylor series approach involves comparing y~Ol
=g~Ol (B), y~ll =gP) (B), . .• , y~n) = g~n) (B), . . .(3)
=
for i 1, .. . , m in order to solve for f). The number of solutions for f) defines the identifiability structure of the system. 2.4 Pohjanpalo 's Rank Test
2.1 General System
For some BO E np , consider the infinite J acobian matrix of derivatives of the Taylor series coefficients of y with respect to the parameter vector e, i.e.
Consider a system of the form x(t)
= f(x , 'U , B, t,)
Y (t)
=g (x,B , t)
(1)
x(O)=xo(B)
where x E !Rn, y E !Rm , u E U (an appropriate space of input functions : usually piecewise continuous functions on an interval, i.e. UfO, T]) and f) E np c !RP where np is the feasible parameter region and f) is a p-dimensional parameter vector. The functions f and 9 are assumed to be analytic in each component.
8g(0) (BO)
og(O) (BO)
oB I og(l) (BO)
oBp og(l) (BO)
oil}
oep
G (eO) =
2.2 Structural Identifiability The identifiability of system (1) is examined at parameter values f) E np in the experiments specified by (xo, UfO , T]), where the initial conditions Xo are well defined once e is selected (Chappell et al. , 1990) . The notation f) - B denotes that the parameter values e,B E np are indistinguishable in the experiments (xo , UfO, T]) . System (1) is defined to be globally identifiable at f) E Dp if f) - B,e E np , implies e B, and it is defined as locally identifiable if there exists an open neighbourhood nn of e in np such that e - G for f) E nn implies f) = B. These definitions can be extended generically so that (1) is globally (locally) structurally identifiable if it is globally (locally) identifiable at all f) E np , except (at most) the points of a subset of zero measure in np. Furthermore (1) is unidentifiable if the system is neither locally or globally identifiable.
=
2.3 Taylor Series Approach
The Taylor series approach (Pohjanpalo, 19i8) examines for the structural identifiability of a
248
(4)
og(k) (eO)
og(k) (eO)
8e 1
oep
Pohjanpalo's result (Pohjanpalo, 1982) states that the parameters e = (f)1 , f)2, .. • ,ep ) .are locally identifiable if and only if G (eO) has rank p. This follows from the Rank theorems provided by Rothenburg (Rothenburg, 1971) . In order to extend the result of Pohjanpalo to consider unidentifiable systems the Implicit Function theorem , and a corollary of the Rank theorem are required (Marsden and Hoffman , 1993) . 3. EXISTENCE AND DETERMINATION OF A REPARAMETERISATION This section presents two theorems. Theorem 1 proves that a reparameterisation exists if the J acobian matrix is rank deficient and Theorem 2 gives conditions which a minimal set of parameters must satisfy such that they can be used to construct a reparameterisation of the original system. Theorem 1. Let Y(f)I, " " Bp , t) be a function whose Taylor series expansion gIVes rise to the coefficients 91 (e), ... , 9p( e) where 9i : A C !RP ---+ !R for i = 1, . . . , p. Define G = (91 , . '" 9p) · Suppose that the Jacobian matrix DG(B) with respect to e has rank q( < p) for all e in a neighbourhood
of 0° E A then the function Y may be locally reparameterised in terms of a set of q of the Taylor series coefficients, . {
=9jl' .. . • <:>q =9J.}
(S)
G(8 1 • .•• , 8p) = G(OI •. . . •
l.e.
(6)
Furthermore the reparameterised system is locally identifiable. Remark : The reparameterisation deriving from the Implicit Function theorem is only guaranteed in a neighbourhood of 0° . As an extension to Theorem 1 let us consider the conditions on 4> for a reparameterisation to exist, i.e what is a suitable (minimal) set of " new" parameters (4)1 ,,,,, 4>q) with which the system can be reparameterised.
5. EXAMPLES To illustrate the methodology presented for reparameterising unidentifiable systems two examples are considered. Example 1 is a model used in Positron Emission Tomography (PET) to characterise ligand binding. The reparameterisation derived in this example puts into a more rigorous context the reparametrisation that has previously been applied empirically (Cunningham et al., 1991). The second example concerns a nonlinear model for microbial growth in a batch reactor. For the purpose of the examples the parameters are assumed to be dimensionless , the dimensions having no effect on the reparameterisation process.
Theorem 2. Given G(81 • .. •• 8p ) = G(
l .... .
(i)
as in Theorem 1 where q < p such that rank (DG(OO)) = q and let the kernel , NG , of (DG(OO)) be given by, NG
= {8 E ~PIDG(8o) . 8 = a}.
(8)
Now given
-
flq C ~q
(9)
suppose that (10 )
where <, > is the Euclidean inner product on and that rank(D qi(60)) = q
~P,
5 .1 Example 1 : Reference Region Model (ll)
A reference region model allows for the compartmental analysis of ligand-receptor interactions without the need for arterial blood sampling (Cunningham et aI., 1991) . The method relies on the fact that there exists a tissue region in which little or no specific binding is present. A region of interest is defined on this reference region such that the kinetics of the free ligand can be obtained. The compartmental model is described by the following set of ordinary linear differential equations ,
then there exists a (; such that G(8 1 . .. .. 8p )
4. METHOD OF
= G(61 ... . ,
q parameters. The final stage involves finding a state space transformation of the original differential equations such that the system may be reparameterised in terms of the new parameter set rP . This step is generally nontrivial and is performed by inspection. By choosing a new parameter set rP involving only simple functions of the original parameters the process of finding a state space transformation is often simplified. A close inspection of the orginal system equations may also prove useful when deriving the most appropriate new parameter set q) . A summary of the reparameterisation process is given in (Gunn and Chappell, 1997).
(12)
REPARA~IETERISATION
The nullspace, N G , of the J acobian matrix G is calculated and it is spanned by, say, (13)
where the vectors nj for i = 1 , ... , P - q are p dimensional. The locally identifiable parameters for the reparameterisation may be calculated from the conditions given in Theorem 2, i.e. the new parameters 4> = (q)I , .. . , q)q) must satisfy, oqii o i T [ < [ -0 , ... , - ] ,nl ... . , n p _ q ] > = 81
08 p
O.
(H)
for i = 1, . . .• q. This produces a partial differential equation whose solution yields the required functions q)i. For the case p - q 1 the partial differential equation is given by
yet)
=
oi 08 1
-nI,
O i + ... + -08 n l . = O. p
= C,(t ) + Co(t)
where C,../(t), C/(t), Cb(t), Ca(t) are the concentrations in the reference region , the free and nonspecific binding pool, the specifically bound pool, and the plasma respectively. J{1 is flow x extraction, k2 is the efflux rate constant, k3 is the product of the association rate constant and the
(IS )
Once q solutions have been found such that rank (D4» =q then by Theorem 2 there exists a reparameterisation of the system 'involving just these
249
maximum specific binding capacity (tracer conditions only), k4 is the dissociation rate constant, J{s is flow x extraction in the reference region and kG is the effiu.x rate constant for the reference region.
dC~.t(t) • d:
;."13k C
1
0
0
01
00
(21)
.
1
It follows that the matrix,G (BD), has rank 5, since using Vajda's result (Vajda, 1982) concerning structural identifiability of linear systems, we know that there are at most 5=(2n-1) independent equations that can be obtained from the Taylor coefficients (where n=3 is the dimension of the state space of the system). Using Theorem 1 there is l=(p- q) redundant parameter and there exists a reparameterisation of the system in terms of just five parameters which is locally identifiable.
(t) re!
1
0
Transformation of the system equations (16) allows us to remove the presence of the term Ca(t). Rearrangement of the first equation gives,
Ca(t) =
°° °° °° °°° -~ °° °°1 ° f°o ° ° l
G' (9°) =
(17)
\5
Substituting this equation into the second of (16) eliminates Ca(t) yielding,
5.1.1. Reparameterisation of the System dCJ(t) KI dCreJ(t) -d-t-;:; Ks dt (k2
+
Klk o ~CreJ(t)-
Calculation of the nullspace of G
+ k3) CJ(t) + k4Cb(t)
(BD)
yields,
[(1 ] n;:; [ K5,0,0,0,1,0 "
(22)
By application of Theorem 2 the locally identifiable parameter groupings can be derived by solving for <{J(B) in equation (10) yielding,
(18)
Kl dtP; Ks dK
1
Cre!(t) may now be treated as the input for this system. The system is examined for a constant input,
dtP;
+ dK5
;:; 0.
(23)
Solution of the differential equation yields the following simple locally identifiable parameter combinations, tPI =
~: ,tP2;:; k2.tP3 = k3,tPi = k4 .rP5 = k6
=
(24)
=
where
This system can be considered as having zero input and non-zero initial conditions, i.e.
dCJ(t) = tPI dCreJ(t) dt
+ tPltPSCrej(t)-
dt (4)2
+ tP3)C j (t) + tPiCb(t) (25)
(20)
Note that the reparameterisation is not unique and that the parameterisation of the system has been reduced from p (= 6) to q (= 5) parameters. The biologically significant parameters are k3 and k4 which characterise the ligand receptor binding and because the reparameterisation allows for the estimation of these parameters it is applicable to experimental data.
The Taylor series of the observation function y(t) expanded about t = 0 can now be calculated. The first five Taylor series coefficients are obtained and an appropriate submatrix of the J acobian matrLx, G (BD), consisting of the partial derivatives of the Taylor series coefficients with respect to the parameter vector B = (J{ 1 , k2, k3, k4' J{5, kG), can then be generated. This process was performed using MATHEMATICA. A subsequent row reduction of this J acobian submatrLx yields,
5.2 Example 2 : Batch Reactor Model
A model of microbial growth in a batch reactor was introduced by Holmberg (Holmberg, 1982) and further analysed by Chappell and Godfrey (Chappell and Godfrey, 1992). The batch reactor model is defined by the equations,
250
= ~m,3 =Kd,
=
(26) 11 (t)
=
ZI
(t)
=0 Z2 (0) = 0 ZI
(0)
where Zt and Z2 are the concentrations of microorganisms and growth limiting substrate respectively, ~m is the maximum velocity of the reaction, K. is the Michaelis-Menten constant, Y is the yield coefficient, Kcl is the decay rate coefficient, bt is the initial concentration of micro-organisms and b2 is the initial concentration of substrate.
(31)
giving the following locally identifiable system ,
The experiment considered is for impulsive inputs to both compartments, Ut (t) = 6 (t) and U2 (t) = c(t) , with observation of Zt . Thus the system can be considered to have zero input to both compartments and non-zero initial conditions, i.e. Ut (t) 0, U2 (t) 0, Zl (0+) b1 , Z2 (0+) = b2 . The first five Taylor series coefficients of the observation, calculated about t = 0, are obtained. Now for some ea E f2p, consider the J acobian matri..x of derivatives of the Tavlor series coefficients of y with respect to the parameter vector e=(~m, K., Kd, bt , Y, b2 ). An appropriate submatrix of the Jacobian matri..x is obtained using MATHEMATICA as the expressions involved are complicated and lengthy and by a process of row reduction this reduces to ,
=
G'
(0°)
=
r
=
=
~ ~ ~ ~ ~ -9<'1
0 0 1 0 0 ~~ 000100 00001
.
(32) 11 (t) ZI
Z2
(27)
(28)
Functions
Y ai
ai
( O.
(0)";: 0
=
=
A method has been provided for reparameterising unidentifiable systems using the Taylor series approach. Theorem 1 presented here states that the rank of the (possibly infinite) J acobian matrLx determines the number of locally identifiable parameters. Theorem 2 gives a sufficient condition for a set of parameters to provide a reparameterisation. A simple method for determining locally identifiable parameter groupings based on the nullspace of the Jacobian matrix has been presented. A reparameterisation of the system is obtained by performing a state space transformation which converts the system into a form in which only identifiable parameter groupings occur. As the method is based on the Taylor series approach it is applicable to linear and nonlinear systems and also to systems with specified inputs (Chappell et al., 1990). In order to calculate the minimal form of the system it is necessary to know the
For the Jacobian matri..."{ G (eO) the nullspace is spanned by the vector,
b;'" aK, - ~ ay + 8~ =
(t)
6. DISCUSSION
5.2.1. RepaTameterisation of the System
y]
ZI
(0) = 0
=
i;
K. = [ 0 , - , 0, 0, - - , 1 . b~ b~
=
where
Hence Rank G (eO) = 5. This model is unidentifiable (Chappell and Godfrey, 1992) and hence the rank of the infinite Jacobian matrix must be less than 6. Using this result and the previous calculations the infinite Jacobian must have rank 5 almost always which implies that q = 5. Using Theorem 1 there is l=(p-q) redundant parameter and there exists a reparameterisation of the system in terms of just 5 parameters which is locally identifiable.
n
=
where
29)
Some possible solutions for
251
rank of the (possibly) infinite Jacobian matrLx (however, for linear models the Jacobian matrLx is finite (Vajda, 1982)). If an infinite Jacobian matrix is rank deficient by 1, this may be deduced by first showing that the system is unidentifiable (e.g. by employing, say, the Similarity transform approach) which by Pohjanpalo's result proves the J acobian must be rank deficient by at least l. Furthermore if a subset of the infinite Jacobian has the necessary rank the J acobian matrLx can be shown to be rank deficient by 1. In cases where the rank of the J acobian is further deficient it may only be possible to get a lower bound on the rank of this matrix. The difficulty arises from the fact that for general nonlinear systems there is no theoretical bound on the number of Taylor series coefficients which need to be considered in order that the J acobian matrix has maximal rank. To reparameterise the system, locally identifiable parameter combinations are established via examination of the nullspace of the J acobian. Then a state space transformation is searched for such that the resulting system's parameters only occur in locally identifiable groupings. It should be stressed that the reparameterisation obtained is not unique, nor is there a specific calculation procedure to follow. Although a state space transformation is usually required to complete the analysis this is not always the case. The construction of a suitable state space transformation is usually non-trivial , but this process can often be greatly simplified by deriving a new (minimal) ~arame ter set which involves only simple functions of the original parameters. Another recommended procedure is to consider the groupings of the parameters in the original system equations when choosing a suitable new (minimal) parameter set. By performing the reparameterisation a robust set of equations is obtained for the experiment considered. The reparameterised model is locally identifiable, by the theory gi ....en , and in addition may be globally identifiable. This method should be applied to models in order to assess whether the biologically significant parameters may be estimated from the experiment considered, i.e. are the biologically significant parameters represented in the reparameterisation ? In practice this method may be implemented easily using any modern symbolic manipulation package. There are two parts of the process which may prove non-trivial; the determination of the rank of the Jacobian and the construction of a suitable state space transformation. Non-trivial unidentifiable systems that have so far been encountered are only rank deficient by one permitting the rank to be determined (see above). The construction of the appropriate state space transformations is normally achieved by inspection and the complexity varies from system to system.
252
7. REFERENCES Chappell, M. J. and K. R. Godfrey (1992). Structural identifiability of the parameters of a nonlinear batch reactor model. M ath Biosci 108, 241-25l. Chappell, M. J., K. R. Godfrey and S. Vajda (1990) . Global identifiability of the parameters of nonlinear systems with specified inputs: a comparison of methods. Math Biosci 102,41-73. Cobelli, C. and G. Toffolo (1987) . Theoretical aspects and practical strategies for the identification of unidentifiable compartmental systems. In: Identifiability of Parametric Models (E. Waiter, Ed.). pp . 85-91. Pergamon Press. Cunningham, V. J., S. P. Hume, G. R. Price, R. G. Ahier, J . E. Cremer and A. K. P. Jones (1991). Compartmental analysis of diprenorphine binding to opiate receptors in the rat in vivo and its comparison with equilibrium data in vitro. J Cereb Blood Flow Metab 11,1-9. Godfrey, K. R. and III DiStefano, J. J. (1987). Identifiability of model parameters. In: Identifiability of Parametric Models (E. Waiter, Ed.). pp. 1-20. Pergamon Press. Gunn, R . N. and M. J . Chappell (1997) . Unidentifiable systems : A methodology for the determination of minimal locally identifiable reparameterisations using the taylor series approach. submitted. Math Biosci. Holmberg, A. (1982). On the practical identifiability of microbial growth models incorporating michaelis-menten type nonlinearities. Math Biosci 62 , 23-43 . Marsden, J. E. and M. J. Hoffman (1993) . Elementary Classical Analysis. second ed .. W. H. Freeman and Company. Milanese, M. and N. Sorrentino (1978) . Decomposition methods for the identifiability analysis oflarge systems. International Journal of Control 1,71-79. Pohjanpalo, H. (1978). System identifiability based on the power series expansion of the solution. Math Biosci 41,21-33. Pohjanpalo, H. (1982). Identifiability of deterministic differential models in state space. Technical report. Technical Research Centre of Finland . Rothenburg, T. J. (1971). Identification in parametric models. Econometrica 39(3), 577-59!. Vajda, S. (1982) . Structural identifiability of linear, bilinear, polynomial and rational systems. In: 9th IFAC World Congress, Budapest, Hungary. Waiter, E. and Y. Lecourtier (1981). Unidentifiable compartmental models: what to do ? Math Biosci 56, 1-25.