This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

11institutetext: Ph.D. Program in Computer Science, CUNY Graduate Center, New York, USA, 11email: [email protected]
22institutetext: Department of Mathematics, CUNY Queens College and Ph.D. Programs in Mathematics and Computer Science, CUNY Graduate Center, New York, USA; 22email: [email protected]
33institutetext: LIX, CNRS, École Polytechnique, Institute Polytechnique de Paris, France
33email: [email protected]

Web-based Structural Identifiability Analyzerthanks: This work was partially supported by the NSF grants CCF-1564132, CCF-1563942, DMS-1760448, DMS-1853650, and DMS-1853482, and by the Paris Ile-de-France Region.

Ilia Ilmer 11    Alexey Ovchinnikov 22    Gleb Pogudin 33
Abstract

Parameter identifiability describes whether, for a given differential model, one can determine parameter values from model equations. Knowing global or local identifiability properties allows construction of better practical experiments to identify parameters from experimental data. In this work, we present a web-based software tool that allows to answer specific identifiability queries. Concretely, our toolbox can determine identifiability of individual parameters of the model and also provide all functions of parameters that are identifiable (also called identifiable combinations) from single or multiple experiments. The program is freely available at https://maple.cloud/app/6509768948056064.

Keywords:
Structural identifiability Identifiability software Differential algebra

1 Introduction and Related Work

A parameter is said to be structurally globally identifiable if, given the input and output of the experiment, one can uniquely recover the parameter’s value in the generic case. If the recovered value is not unique but comes from a finite collection, then we say that such a parameter is locally identifiable. Otherwise, the parameter is called non-identifiable. In the latter case, one wonders if there is a function of that parameter that is identifiable. This is useful in several ways, for instance, it can mitigate the issue of non-identifiability of some parameters [13].

There is a variety of installable packages that deal with parameter identifiability, see, for instance [1, 4, 10, 17, 21]. For a more detailed overview of these, see [5, 8] and references therein. A general overview of solving parameter identifiability problems was presented, for instance, in [12, 14, 20]. Among the available identifiability software, SIAN [7] written in Maple111For a Julia implementation, see https://github.com/alexeyovchinnikov/SIAN-Julia is typically the fastest one for assessing global identifiability of individual parameters (see, e.g., [7, Table 1]). In the case of the lack of identifiability, one may want to find which functions of the parameters are identifiable. For this task, DAISY software [1] implemented in Reduce can be used (under some assumptions, see [15, Remark 3] and [13]). This state of affairs may be inconvenient for the user because

  1. 1)

    the features of interest are scattered among different packages;

  2. 2)

    packages may require proprietary (Maple) or less popular (Reduce) software and may not be available for commonly used OS (DAISY is not available for the UNIX-type systems);

  3. 3)

    finally, the packages should be installed.

These issues have been partially addressed by a web-based tool called COMBOS [11] (and its recent refinement COMBOS 2 for linear systems [9]). However, the backend algorithm appears to be less efficient than SIAN [7, Table 1], and it relies on the same assumption on the input model as DAISY.

Our main contribution is a web-based toolbox hosted on Maple Cloud for assessing structural identifiability built upon SIAN and recent software for computing identifiable functions of parameters [13] which uses the Boulier’s BLAD software package [2] incorporated into the Maple’s Differential Algebra package. The key features are

  1. 1)

    efficiency. We use SIAN for assessing identifiability of individual parameters efficiently. For computing identifiable functions, we use the code from [13] which we speed up by exploiting the results of the computation performed by SIAN.

  2. 2)

    versatility. The toolbox allows assessing local and global identifiability of the parameters and initial conditions and compute the identifiable functions in parameter both in the single- and multi-experiment setup. We do not make any assumptions on the input system unlike DAISY or COMBOS.

  3. 3)

    availability. The toolbox is a web app, so it can be used in a browser in one click and does not require installing anything.

In Section 3, we outline several scenarios in which our application is essential for assessing identifiability of parameters and parameter combinations. We also illustrate the speedup achievable using output from each of its parts. The web-application222The Maple implementations of each underlying algorithm are available on GitHub at https://github.com/pogudingleb/SIAN and https://github.com/pogudingleb/AllIdentifiableFunctions. can be used at https://maple.cloud/app/6509768948056064 and is also available for download.

2 Input-output specification

Let us define the specific form of state-space input ODE that our application accepts.

Definition 1 (Model in the state-space form)

A model in the state-space form accepted by the application is a system

𝚺:={𝐱=𝐟(𝐱,𝝁,𝐮),𝐲=𝐠(𝐱,𝝁,𝐮),𝐱(0)=𝐱,\mathbf{\Sigma}:=\begin{cases}\mathbf{x}^{\prime}&=\mathbf{f}(\mathbf{x},\boldsymbol{\mu},\mathbf{u}),\\ \mathbf{y}&=\mathbf{g}(\mathbf{x},\boldsymbol{\mu},\mathbf{u}),\\ \mathbf{x}(0)&=\mathbf{x}^{\ast},\end{cases}

where 𝐟=(f1,,fn)\mathbf{f}=(f_{1},\dots,f_{n}) and 𝐠=(g1,,gn)\mathbf{g}=(g_{1},\dots,g_{n}) with fi=fi(𝐱,𝝁,𝐮)f_{i}=f_{i}(\mathbf{x},\boldsymbol{\mu},\mathbf{u}), gi=gi(𝐱,𝝁,𝐮)g_{i}=g_{i}(\mathbf{x},\boldsymbol{\mu},\mathbf{u}) are rational functions over the field of complex numbers \mathbb{C}.

The vector 𝐱=(x1,,xn)\mathbf{x}=(x_{1},\dots,x_{n}) represents the time-dependent state variables and 𝐱\mathbf{x}^{\prime} represents the derivative. The vector-function 𝐮=(u1,,us)\mathbf{u}=(u_{1},\dots,u_{s}) represents the input variable. The mm-vector 𝐲=(y1,,yn)\mathbf{y}=(y_{1},\dots,y_{n}) represents the output variables. The vector 𝝁=(μ1,,μλ)\boldsymbol{\mu}=(\mu_{1},\dots,\mu_{\lambda}) represents the parameters and 𝐱=(x1,,xn)\mathbf{x}^{\ast}=(x_{1}^{\ast},\dots,x_{n}^{\ast}) defines initial conditions of the model.

Below we specify the input format and possible outputs of our toolbox. Note that while used in descriptions below, some outputs, such as number of solutions for each parameter, are not listed here for brevity. The app also provides additional logs for debugging purposes. In Appendix 0.B, we provide more specification examples.

In:A model in state-space form, see Definition 1.\displaystyle\begin{array}[]{rl}\textbf{In:}&\text{A model in state-space form, see \lx@cref{creftypecap~refnum}{def1}.}\end{array}
Out:Globally:Globally identifiable parameters, that isones uniquely recoverable for a given system.Locally not Globally:Locally but not globally identifiableparameters, with finitely many recoverable values.Non-Identifiable:Non-identifiable parameters, these can haveinfinitely many values.Single-Experiment:Single-Experiment identifiable functions ofparameters, i.e. identifiable from k1 experiments.Multi-Experiment:Multi-Experiment identifiable functions ofparameters, i.e. identifiable from kβ experiments.β:Bound on the number of experiments.\displaystyle\begin{array}[]{rll}\\ \textbf{Out:}&\textit{Globally}:&\text{Globally identifiable parameters, that is}\\ &&\text{ones uniquely recoverable for a given system.}\\ &\textit{Locally not Globally}:&\text{Locally but not globally identifiable}\\ &&\text{parameters, with finitely many recoverable values.}\\ &\textit{Non-Identifiable}:&\text{Non-identifiable parameters, these can have}\\ &&\text{infinitely many values.}\\ &\textit{Single-Experiment}:&\text{Single-Experiment identifiable functions of}\\ &&\text{parameters, i.e. identifiable from $k\leq 1$ experiments.}\\ &\textit{Multi-Experiment}:&\text{Multi-Experiment identifiable functions of}\\ &&\text{parameters, i.e.\ identifiable from $k\leq\beta$ experiments.}\\ &\beta:&\text{Bound on the number of experiments.}\end{array}

Note that the single- and multi-experiment identifiable combinations returned by the app generate all single- and multi-experiment functions of parameters, respectively. We return them in the algebraically simplified form. In addition, the app reports number of solutions per each globally or locally identifiable parameter, which is not explicitly reflected here due to space limitations.

3 Use Cases for Structural Identifiability Toolbox

3.1 Globally Identifiable Example (two-species competition model)

Let us consider a simple two-species competition model based with logistic growth in homogeneous environment and assume that we are interested in identifiability properties of all parameters and initial conditions:

{x1=r1x1(1x1+x2k1),x2=r2x2(1x1+x2k2),y1=x1,y2=x2\begin{cases}x_{1}^{\prime}=r_{1}x_{1}\left(1-\tfrac{x_{1}+x_{2}}{k_{1}}\right),\\ x_{2}^{\prime}=r_{2}x_{2}\left(1-\tfrac{x_{1}+x_{2}}{k_{2}}\right),\\ y_{1}=x_{1},~{}~{}y_{2}=x_{2}\end{cases}

with population densities x1,x2x_{1},x_{2} being time-dependent state variables, and intrinsic growth rates r1,r2r_{1},r_{2} and carrying capacities k1,k2k_{1},k_{2} being constant. To run the toolbox for this system, we would write the following into the input field:

In: diff(x1(t),t) = r1*x1(t)*(1 - (x1(t) + x2(t))/k1);diff(x2(t),t) = r2*x2(t)*(1 - (x1(t) + x2(t))/k2);y1(t) = x1(t);y2(t) = x2(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}~{}}&\texttt{diff(x1(t),t) = r1*x1(t)*(1 - (x1(t) + x2(t))/k1);}\\ &\texttt{diff(x2(t),t) = r2*x2(t)*(1 - (x1(t) + x2(t))/k2);}\\ &\texttt{y1(t) = x1(t);}\\ &\texttt{y2(t) = x2(t)}\end{array}
Out: Globally:[x1(0), x2(0), r1, r2, k1, k2]Locally not Globally:[]Non-Identifiable:[]\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[x1(0), x2(0), r1, r2, k1, k2]}\\ &\textit{Locally not Globally}:&\texttt{[]}\\ &\textit{Non-Identifiable}:&\texttt{[]}\end{array}

To determine the identifiability for this model, we keep default “Check global/local identifiability” and “Print Number of Solutions” options on. After entering the system and running the application, the output field contains the results. In this model, all parameters and initial conditions are globally (and locally) identifiable. One can now proceed to data collection and further experiments.

3.2 Locally Identifiable Model (SIRS model with forcing)

Consider an example of a seasonal epidemic model with a periodic forcing term:

{s=μμsb0(1+b1x1)is+gr,i=b0(1+b1x1)is(ν+μ)i,r=νi(μ+g)r,x1=Mx2,x2=Mx1,y1=i,y2=r.\begin{cases}s^{\prime}=\mu-\mu s-b_{0}(1+b_{1}x_{1})i\cdot s+g\cdot r,\\ i^{\prime}=b_{0}(1+b_{1}x_{1})i\cdot s-(\nu+\mu)i,\\ r^{\prime}=\nu i-(\mu+g)r,\\ x_{1}^{\prime}=-Mx_{2},\\ x_{2}^{\prime}=Mx_{1},\\ y_{1}=i,\ \ y_{2}=r.\end{cases}

The model is taken from [3] and is built into the application as one of the illustrating examples. Assume that we are interested in identifiability of parameters of this model. Without changing default settings, running the application yields the result of b1,x1(0),x2(0)b_{1},x_{1}(0),x_{2}(0) being unidentifiable, and b0,g,μ,ν,s(0),i(0),r(0)b_{0},g,\mu,\nu,s(0),i(0),r(0) as globally identifiable. At the same time, we observe that MM which defines oscillation of the term x1x_{1} is the only parameter identifiable locally, not globally. By checking the number of solutions, we see that only two can be found for MM with probability p=0.99p=0.99. Since MM represents the oscillation frequency, it is assumed to be positive in practice, hence globally identifiable. Note that we only needed a single section of the app and the result has been obtained in about 7.2 seconds.

3.3 Identifiable Combination of Non-Identifiable Parameters (tumor targeting)

In this example, we consider system 3 from [16, Section 3] with unknown initial conditions. The example describes a compartmental model describing tumor targeting with antibodies, see [18]. To arrive at the system below, we suppose equations (B) and (D) are identically zero and that 5V36V3=1\frac{5V36}{V3}=1. The functions xi,i=1,,5x_{i},i=1,\dots,5 represent concentrations, ki,i=3,,7k_{i},i=3,\dots,7 and a,b,da,b,d represent rate constants.

{x1=(k3+k7)x1+k4x2,x2=k3x1(k4+(a+bd)k5)x2+k6(x3+x4)+k5x2(x3+x4),x3=ak5x2k6x3k5x2x3,x4=bdk5x2k6x4k5x2x4,x5=k7x1,y1=x5.\begin{cases}x_{1}^{\prime}=-(k_{3}+k_{7})x_{1}+k_{4}x_{2},\\ x_{2}^{\prime}=k_{3}x_{1}-(k_{4}+(a+bd)k_{5})x_{2}+k_{6}(x_{3}+x_{4})+k_{5}x_{2}(x_{3}+x_{4}),\\ x_{3}^{\prime}=ak_{5}x_{2}-k_{6}x_{3}-k_{5}x_{2}x_{3},\\ x_{4}^{\prime}=bdk_{5}x_{2}-k_{6}x_{4}-k_{5}x_{2}x_{4},\\ x_{5}^{\prime}=k_{7}x_{1},\\ y_{1}=x_{5}.\end{cases}

For this model, after computing identifiability properties using SIAN, we observe that everything except parameters a,b,d,x3(0),x4(0)a,b,d,x_{3}(0),x_{4}(0) is globally or locally identifiable. To investigate further, we consider computation with “Compute Identifiable Combinations” option turned on. Running the program with this additional setting, we see that while parameters a,b,da,b,d are not identifiable, their combination a+bda+bd can be identified from at most one experiment. This is especially beneficial since one can connect the meaning of expression a+bda+bd to the overall biological sense of the model’s underlying phenomenon. For instance, in the original paper [18], constant aa and a product bdbd may be attributed to total binding sites on normal tissue and number of binding sites on tumor making a+bda+bd the total number of binding sites in the system. Further, one could apply a substitution of the form x^3=x3+x4,p^=a+bd\widehat{x}_{3}=x_{3}+x_{4},\widehat{p}=a+bd so that in the new system we only have equations for x1,x2,x^3,x5x_{1},x_{2},\widehat{x}_{3},x_{5} and the parameter combination a+bda+bd will now be globally identifiable as a parameter p^\widehat{p}.

3.4 System with a Non-Identifiable Parameter (Lotka-Volterra model)

Let us consider the following Lotka-Volterra model

{x1=ax1bx1x2x2=cx2+dx1x2y=x1\begin{cases}x_{1}^{\prime}=ax_{1}-bx_{1}x_{2}\\ x_{2}^{\prime}=-cx_{2}+dx_{1}x_{2}\\ y=x_{1}\end{cases} (1)

By running the application for (1) using only SIAN, we see that parameter bb and initial condition x2(0)x_{2}(0) are non-identifiable, and the parameters a,b,da,b,d and the initial condition x1(0)x_{1}(0) are globally identifiable. Furthermore, since aa is identifiable and x1x_{1} is observed, from the first equation we conclude that bx2(0)=ax1(0)/x1(0)bx_{2}(0)=a-x_{1}^{\prime}(0)/x_{1}(0) is identifiable. This implies that we have an output-preserving scaling transformation bλbb\to\lambda b, x2x2/λx_{2}\to x_{2}/\lambda. Therefore, the reparametrization x^2:=bx2\hat{x}_{2}:=bx_{2} makes the model globally identifiable.

3.5 Refining Multi-Experiment Identifiability Bound
(slow-fast ambiguity in a chemical reaction network)

Consider the following system:

{xA=k1xA,xB=k1xAk2xB,xC=k2xB,eA=eC=0,y1=eAxA+eBxB+eCxC,y2=xC,y3=eA,y4=eC.\begin{cases}x_{A}^{\prime}=-k_{1}x_{A},\\ x_{B}^{\prime}=k_{1}x_{A}-k_{2}x_{B},\\ x_{C}^{\prime}=k_{2}x_{B},\\ e_{A}^{\prime}=e_{C}^{\prime}=0,\\ y_{1}=e_{A}x_{A}+e_{B}x_{B}+e_{C}x_{C},\\ y_{2}=x_{C},~{}y_{3}=e_{A},~{}y_{4}=e_{C}.\end{cases}

This model is based on a kinetic reaction Ak1Bk2CA\xrightarrow[]{k_{1}}B\xrightarrow[]{k_{2}}C from [19] and has an extra output equation y2y_{2}. The functions xA,xB,xCx_{A},~{}x_{B},~{}x_{C} are concentrations and eA,eBe_{A},~{}e_{B} with constant eCe_{C} represent molar extinction coefficients. In addition, parameters include unknown rate coefficients k1,k2k_{1},~{}k_{2}. The application reports global identifiability for xC(0),eA(0),eC(0)x_{C}(0),~{}e_{A}(0),~{}e_{C}(0) and local identifiability for everything else.

It is then of interest to check identifiable parameter combinations. The app reports single-experiment identifiability for k1k2,k1+k2k_{1}k_{2},~{}k_{1}+k_{2}. This implies that the parameters k1k_{1} and k2k_{2} are identifiable up to a permutation, so it is possible to infer the reaction rates from an experiment but not which rate corresponds to which reaction. Interestingly, the app reports that eB,k1,k2e_{B},~{}k_{1},~{}k_{2} become globally identifiable if one performs at most 3 experiments. Can we do better? To answer this, we turn on the option “Try to Refine Bound” with default number of refining attempts being 4. As a result, the app reports a new bound for the number of experiments being 2.

Let us illustrate this point in another way. Recall that we can tell SIAN to consider multiple copies of the system when analyzing identifiability. In this mode, SIAN does not output initial conditions for brevity. We observed that the refined bound for parameters eB,k1,k2e_{B},k_{1},k_{2} was 2. If we set the “Number of experiments (copies of the input system)” to 2, SIAN yields global identifiability of eB,k1,k2e_{B},k_{1},k_{2}, which verifies our earlier finding. Moreover, turning off “Attempt Bypass using SIAN” option in the search for combinations, we observe that the application still returns eB,k1,k2e_{B},k_{1},k_{2} as identifiable with 3 experiments, however, single experiment check overwrites this result, yielding bound of 1.

Acknowledgements

We are grateful to Joseph DiStefano III, Jürgen Gerhard, John May, Maria Pia Saccomani, and Eduardo Sontag for fruitful discussions, useful feedback, and technical assistance.

References

  • [1] Giuseppina Bellu, Maria Pia Saccomani, Stefania Audoly and Leontina D’Angiò “DAISY: A new software tool to test global identifiability of biological and physiological systems” In Computer methods and programs in biomedicine 88.1 Elsevier, 2007, pp. 52–61 URL: https://doi.org/10.1016/j.cmpb.2007.07.002
  • [2] François Boulier and Équipe Calcul Formel “BLAD–Bibliothèques Lilloises d’Algèbre Différentielle”, 2014 URL: https://cristal.univ-lille.fr/~boulier/pmwiki/pmwiki.php/Main/BLAD
  • [3] Marcos A Capistrán, Miguel A Moreles and Bruno Lara “Parameter estimation of some epidemic models. The case of recurrent epidemics caused by respiratory syncytial virus” In Bulletin of mathematical biology 71.8 Springer, 2009, pp. 1890 URL: https://doi.org/10.1007/s11538-009-9429-3
  • [4] Oana-Teodora Chiş, Julio R Banga and Eva Balsa-Canto “GenSSI: a software toolbox for structural identifiability analysis of biological models” In Bioinformatics 27.18 Oxford University Press, 2011, pp. 2610–2611 URL: https://doi.org/10.1093/bioinformatics/btr431
  • [5] Oana-Teodora Chiş, Julio R. Banga and Eva Balsa-Canto “Structural Identifiability of Systems Biology Models: A Critical Comparison of Methods” In PLoS ONE 6.11, 2011, pp. e27755 URL: https://doi.org/10.1371/journal.pone.0027755
  • [6] Carsten Conradi and Anne Shiu “Dynamics of Posttranslational Modification Systems: Recent Progress and Future Directions” In Biophysical Journal 114.3, 2018, pp. 507–515 DOI: https://doi.org/10.1016/j.bpj.2017.11.3787
  • [7] Hoon Hong, Alexey Ovchinnikov, Gleb Pogudin and Chee Yap “SIAN: software for structural identifiability analysis of ODE models” In Bioinformatics 35.16, 2019, pp. 2873–2874 URL: https://doi.org/10.1093/bioinformatics/bty1069
  • [8] Hoon Hong, Alexey Ovchinnikov, Gleb Pogudin and Chee Yap “Global identifiability of differential models” In Communications on Pure and Applied Mathematics 73.9 Wiley Online Library, 2020, pp. 1831–1879 URL: https://doi.org/10.1002/cpa.21921
  • [9] Ali Kalami Yazdi, Mehdi Nadjafikhah and Joseph Distefano III “COMBOS2: an algorithm to the input–output equations of dynamic biosystems via Gaussian elimination” In Journal of Taibah University for Science 14.1 Taylor & Francis, 2020, pp. 896–907 URL: https://doi.org/10.1080/16583655.2020.1776466
  • [10] Thomas S Ligon et al. “GenSSI 2.0: multi-experiment structural identifiability analysis of SBML models” In Bioinformatics 34.8 Oxford University Press, 2018, pp. 1421–1423 URL: https://doi.org/10.1093/bioinformatics/btx735
  • [11] Nicolette Meshkat, Christine Er-zhen Kuo and Joseph DiStefano III “On finding and using identifiable parameter combinations in nonlinear dynamic systems biology models and COMBOS: a novel web implementation” In PLoS One 9.10 Public Library of Science, 2014, pp. e110261 URL: https://doi.org/10.1371/journal.pone.0110261
  • [12] Hongyu Miao, Xiaohua Xia, Alan S Perelson and Hulin Wu “On identifiability of nonlinear ODE models and applications in viral dynamics” In SIAM review 53.1 SIAM, 2011, pp. 3–39 URL: https://doi.org/10.1137/090757009
  • [13] Alexey Ovchinnikov, Anand Pillay, Gleb Pogudin and Thomas Scanlon “Computing all identifiable functions for ODE models” In arXiv preprint arXiv:2004.07774, 2020 URL: https://arxiv.org/abs/2004.07774
  • [14] Andreas Raue et al. “Comparison of approaches for parameter identifiability analysis of biological systems” In Bioinformatics 30.10 Oxford University Press, 2014, pp. 1440–1448 URL: https://doi.org/10.1093/bioinformatics/btu006
  • [15] M. Saccomani, S. Audoly and L. D’Angiò “Parameter identifiability of nonlinear systems: the role of initial conditions” In Automatica 39, 2003, pp. 619–632 URL: https://doi.org/10.1016/S0005-1098(02)00302-3
  • [16] Maria Pia Saccomani, Stefania Audoly, Giuseppina Bellu and Leontina D’Angiò “Examples of testing global identifiability of biological and biomedical models with the DAISY software” In Computers in Biology and Medicine 40.4 Elsevier, 2010, pp. 402–407 URL: https://doi.org/10.1016/j.compbiomed.2010.02.004
  • [17] Maria Pia Saccomani, Giuseppina Bellu, Stefania Audoly and Leontina D’Angiò “A new version of DAISY to test structural identifiability of biological models” In International conference on computational methods in systems biology, 2019, pp. 329–334 Springer URL: https://doi.org/10.1007/978-3-030-31304-3˙21
  • [18] Gillian D Thomas et al. “Effect of dose, molecular size, affinity, and protein binding on tumor uptake of antibody or ligand: a biomathematical model” In Cancer research 49.12 AACR, 1989, pp. 3290–3296 URL: http://www.ncbi.nlm.nih.gov/pubmed/2720683
  • [19] S Vajda and H Rabitz “Identifiability and distinguishability of first-order reaction systems” In The Journal of Physical Chemistry 92.3 ACS Publications, 1988, pp. 701–707 URL: https://doi.org/10.1021/j100314a024
  • [20] Alejandro F Villaverde “Observability and structural identifiability of nonlinear biological systems” In Complexity 2019 Hindawi, 2019 URL: https://doi.org/10.1155/2019/8497093
  • [21] Alejandro F Villaverde, Antonio Barreiro and Antonis Papachristodoulou “Structural identifiability of dynamic systems biology models” In PLoS computational biology 12.10 Public Library of Science San Francisco, CA USA, 2016, pp. e1005153 URL: https://doi.org/10.1371/journal.pcbi.1005153

Appendix 0.A Details on the underlying algorithms

The application solves two problems: identifiability properties of individual parameters and that of combinations (functions) of parameters. Note that we return generators of the field of all identifiable functions.

The input for both problems follows the same structure where we pass a collection of ODEs and output functions. For querying identifiability of individual parameters and initial conditions we use SIAN [7]. In short, it expresses the Taylor coefficients of output functions in terms of the initial conditions and parameters and checks whether the parameters or initial conditions of interest can be expressed via these coefficients. For better efficiency, this is checked for a randomly sampled solution of the system. The probability that such a solution will exhibit the generic behavior is quantified in [8]. Therefore, the overall algorithm is randomized Monte Carlo, that is the result is guaranteed correct with user-specified probability pp.

To answer the question on identifiability of parameter functions we take advantage of work [13]. Note that our application distinguishes single- and multi-experiment identifiable combinations as opposed to existing methods for identifiable combination queries. The latter is equivalent to having multiple copies of original ODE system sharing the parameters, outputs, and inputs. We also provide a bound on the number of experiments which can be refined by changing ordering of variables in the underlying algorithm.

We compute the input-output equations, that is, differential equations relating inputs, outputs, and parameters of the differential model. Identifiable functions of parameters are then extracted from the coefficients of these equations using methods of differential algebra and computational algebraic geometry, including Gröbner basis computation. To minimize the computational overhead, we take advantage of the Gröbner walk procedure, by changing the order from total degree reverse to pure lexicographic. This algorithm is deterministic or a Monte Carlo probabilistic, depending on how/which of the Gröbner basis implementation is used.

To achieve maximal speed of computation without compromising the functionality of the application, we take advantage of the fact that SIAN is typically faster than the algorithm from [13] and its output can be sometimes used to obtain the output of [13] without further computation. More precisely, if all parameters are reported as globally identifiable with probability pp, then, with the same probability, we report these parameters as their own identifiable combinations and an example of this is presented in Section 0.B.5.

With the current implementation, the application does not support specifying initial conditions, however this functionality is planned for future versions.

Refer to caption
Figure 1: Main view of the application. The dial on the right side indicates whether an app is running. The arrows are clickable and show additional settings for each section of the program as well as documentation.
Refer to caption
Figure 2: Output fields of the application. We present individual parameters’ results and identifiable combinations with the bound separately. The white field on the left displays number of solutions per identifiable parameter.

Appendix 0.B Systems in Structural Identifiability Toolbox Input Form

Below we present input and output form for examples discussed in this paper. The input is shown in the Maple syntax form. The toolbox also supports a different input format:

dx1/dt = a*x1 + x2*b + u(t);
dx2/dt = x2*c + x1;
y = x2

where inputs u(t) are required to have argument written explicitly.

0.B.1 Example from Section 3.2

In: diff(s(t), t) = mu - mu*s(t) - b0*(1 + b1*x1(t))*i(t)*s(t)+ g(t)*r(t);diff(i(t), t) = b0*(1 + b1*x1(t))*i(t)*s(t) - (nu+mu)*i(t);diff(r(t), t) = nu*i(t) - (mu + g)*r(t);diff(x1(t), d) = -M*x2(t);diff(x2(t), d) = M*x1(t);y1(t) = i(t);y2(t) = r(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}~{}}&\texttt{diff(s(t), t) = mu - mu*s(t) - b0*(1 + b1*x1(t))*i(t)*s(t)}\\ &\texttt{+ g(t)*r(t);}\\ &\texttt{diff(i(t), t) = b0*(1 + b1*x1(t))*i(t)*s(t) - (nu+mu)*i(t);}\\ &\texttt{diff(r(t), t) = nu*i(t) - (mu + g)*r(t);}\\ &\texttt{diff(x1(t), d) = -M*x2(t);}\\ &\texttt{diff(x2(t), d) = M*x1(t);}\\ &\texttt{y1(t) = i(t);}\\ &\texttt{y2(t) = r(t)}\\ \end{array}
Out: Globally:[b0, g, mu, nu, s(0), i(0), r(0)]Locally not Globally:[M]Non-Identifiable:[b1, x1(0), x2(0)]\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[b0, g, mu, nu, s(0), i(0), r(0)]}\\ &\textit{Locally not Globally}:&\texttt{[M]}\\ &\textit{Non-Identifiable}:&\texttt{[b1, x1(0), x2(0)]}\end{array}

0.B.2 Example from Section 3.3

In: diff(x1(t),t) = -(k3 + k7)*x1(t) + k4*x2(t);diff(x2(t), t) = k3*x1(t) - (k4 + (a + b*d)*k5)*x2(t)+ k6*(x3(t) + x4(t)) + k5*x2(t)*(x3(t) + x4(t));diff(x3(t), t) = a*k5*x2(t) - k6*x3(t) - k5*x2(t)*x3(t);diff(x4(t), t) = b*d*k5*x2(t) - k6*x4(t) - k5*x2(t)*x4(t);diff(x5(t), t) = k7*x1(t);y1(t) = x5(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}~{}}&\texttt{diff(x1(t),t) = -(k3 + k7)*x1(t) + k4*x2(t);}\\ &\texttt{diff(x2(t), t) = k3*x1(t) - (k4 + (a + b*d)*k5)*x2(t)}\\ &\texttt{+ k6*(x3(t) + x4(t)) + k5*x2(t)*(x3(t) + x4(t));}\\ &\texttt{diff(x3(t), t) = a*k5*x2(t) - k6*x3(t) - k5*x2(t)*x3(t);}\\ &\texttt{diff(x4(t), t) = b*d*k5*x2(t) - k6*x4(t) - k5*x2(t)*x4(t);}\\ &\texttt{diff(x5(t), t) = k7*x1(t);}\\ &\texttt{y1(t) = x5(t)}\\ \end{array}
Out: Globally:[k3, k4, k5, k6, k7, x1(0), x2(0), x5(0)]Locally not Globally:[]Non-Identifiable:[a, b, d, x3(0), x4(0)]Single-Experiment:[k3, k4, k6, k7, k5/k7, bd+a]Multi-Experiment:[k3, k4, k6, k7, k5/k7, bd+a]β=1\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[k3, k4, k5, k6, k7, x1(0), x2(0), x5(0)]}\\ &\textit{Locally not Globally}:&\texttt{[]}\\ &\textit{Non-Identifiable}:&\texttt{[a, b, d, x3(0), x4(0)]}\\ &\textit{Single-Experiment}:&\texttt{[k3, k4, k6, k7, k5/k7, bd+a]}\\ &\textit{Multi-Experiment}:&\texttt{[k3, k4, k6, k7, k5/k7, bd+a]}\\ &\beta=&1\end{array}

0.B.3 Example from Section 3.4

In: diff(x1(t), t) = a*x1(t) - b*x1(t)*x2(t);diff(x2(t), t) = -c*x2(t) + d*x1(t)*x2(t);y(t) = x1(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}~{}}&\texttt{diff(x1(t), t) = a*x1(t) - b*x1(t)*x2(t);}\\ &\texttt{diff(x2(t), t) = -c*x2(t) + d*x1(t)*x2(t);}\\ &\texttt{y(t) = x1(t)}\\ \end{array}
Out: Globally:[a, c, d, x1(0)]Locally not Globally:[]Non-Identifiable:[b, x2(0)]Single-Experiment:[a, c, d]Multi-Experiment:[a, c, d]β=1\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[a, c, d, x1(0)]}\\ &\textit{Locally not Globally}:&\texttt{[]}\\ &\textit{Non-Identifiable}:&\texttt{[b, x2(0)]}\\ &\textit{Single-Experiment}:&\texttt{[a, c, d]}\\ &\textit{Multi-Experiment}:&\texttt{[a, c, d]}\\ &\beta=&1\end{array}

0.B.4 Example from Section 3.5

In: diff(xA(t), t) = -k1*xA(t);diff(xB(t), t) = k1*xA(t) - k2*xB(t);diff(xC(t), t) = k2*xB(t);diff(eA(t), t) = 0;diff(eC(t), t) = 0;y1(t) = eA(t)*xA(t) + eB*xB(t) + eC(t)*xC(t);y2(t) = xC(t);y3(t) = eA(t);y4(t) = eC(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}~{}}&\texttt{diff(xA(t), t) = -k1*xA(t);}\\ &\texttt{diff(xB(t), t) = k1*xA(t) - k2*xB(t);}\\ &\texttt{diff(xC(t), t) = k2*xB(t);}\\ &\texttt{diff(eA(t), t) = 0;}\\ &\texttt{diff(eC(t), t) = 0;}\\ &\texttt{y1(t) = eA(t)*xA(t) + eB*xB(t) + eC(t)*xC(t);}\\ &\texttt{y2(t) = xC(t);}\\ &\texttt{y3(t) = eA(t);}\\ &\texttt{y4(t) = eC(t)}\\ \end{array}
Out: Globally:[xC(0), eA(0), eC(0)]Locally not Globally:[eB, k1, k2, xA(0), xB(0)]Non-Identifiable:[]Single-Experiment=[k1k2, k1+k2]Multi-Experiment=[eB, k1, k2]β=3\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[xC(0), eA(0), eC(0)]}\\ &\textit{Locally not Globally}:&\texttt{[eB, k1, k2, xA(0), xB(0)]}\\ &\textit{Non-Identifiable}:&\texttt{[]}\\ &\textit{Single-Experiment}=&\texttt{[k1k2, k1+k2]}\\ &\textit{Multi-Experiment}=&\texttt{[eB, k1, k2]}\\ &\beta=&3\end{array}

0.B.5 Example of speedup with Bypasses

Consider the system from [6]

{x1=k1x1x2+k2x4+k4x6,x2=k1x1x2+k2x4+k3x4,x3=k3x4+k5x6k6x3x5,x4=k1x1x2k2x4k3x4,x5=k4x6+k5x6k6x3x5,x6=k4x6k5x6+k6x3x5,y1=x3,y2=x2.\begin{cases}x_{1}^{\prime}=-k_{1}x_{1}x_{2}+k_{2}x_{4}+k_{4}x_{6},\\ x_{2}^{\prime}=k_{1}x_{1}x_{2}+k_{2}x_{4}+k_{3}x_{4},\\ x_{3}^{\prime}=k_{3}x_{4}+k_{5}x_{6}-k_{6}x_{3}x_{5},\\ x_{4}^{\prime}=k_{1}x_{1}x_{2}-k_{2}x_{4}-k_{3}x_{4},\\ x_{5}^{\prime}=k_{4}x_{6}+k_{5}x_{6}-k_{6}x_{3}x_{5},\\ x_{6}^{\prime}=-k_{4}x_{6}-k_{5}x_{6}+k_{6}x_{3}x_{5},\\ y_{1}=x_{3},\\ y2=x_{2}.\end{cases}

This is an example of a mixed-mechanism network, where the state functions xi(t),i=1,,6x_{i}(t),i=1,\dots,6 are concentrations and the parameters ki,i=1,,6k_{i},i=1,\dots,6 are rate constants. The app returns global and local identifiability for all parameter in under 4 seconds. This is used to conclude that multi-experiment identifiable combinations with the bound of 1 are parameters themselves. If we turn off the “Attempt Bypass” function, the multi-experiment identifiable combinations k1,k3,k5,k6,k2k4+k3k5k2+k3,k2k3k_{1},k_{3},k_{5},k_{6},\frac{-k_{2}k_{4}+k_{3}k_{5}}{k_{2}+k_{3}},k_{2}-k_{3} with bound 1 are returned in 433 seconds. The input form for this example is presented below

In: diff(x1(t), t) =-k1*x1(t)*x2(t) + k2*x4(t) + k4*x6(t);diff(x2(t), t) = k1*x1(t)*x2(t) + k2*x4(t) + k3*x4(t);diff(x3(t), t) = k3*x4(t) + k5*x6(t) - k6*x3(t)*x5(t);diff(x4(t), t) = k1*x1(t)*x2(t) - k2*x4(t) - k3*x4(t);diff(x5(t), t) = k4*x6(t) + k5*x6(t) - k6*x3(t)*x5(t);diff(x6(t), t) =-k4*x6(t) - k5*x6 (t)+ k6*x3(t)*x5(t);y1(t) = x3(t);y2(t) = x2(t)\displaystyle\begin{array}[]{r@{}l}\textbf{In:~{}}&\texttt{diff(x1(t), t) =-k1*x1(t)*x2(t) + k2*x4(t) + k4*x6(t);}\\ &\texttt{diff(x2(t), t) = k1*x1(t)*x2(t) + k2*x4(t) + k3*x4(t);}\\ &\texttt{diff(x3(t), t) = k3*x4(t) + k5*x6(t) - k6*x3(t)*x5(t);}\\ &\texttt{diff(x4(t), t) = k1*x1(t)*x2(t) - k2*x4(t) - k3*x4(t);}\\ &\texttt{diff(x5(t), t) = k4*x6(t) + k5*x6(t) - k6*x3(t)*x5(t);}\\ &\texttt{diff(x6(t), t) =-k4*x6(t) - k5*x6 (t)+ k6*x3(t)*x5(t);}\\ &\texttt{y1(t) = x3(t);}\\ &\texttt{y2(t) = x2(t)}\\ \end{array}
Out: Globally:[k1, k2, k3, k4, k5, k6,x1(0), x2(0), x3(0),x4(0), x5(0), x6(0)]Locally not Globally:[]Non-Identifiable:[]Single-Experiment:[k1, k2, k3, k4, k5, k6]Multi-Experiment:[k1, k2, k3, k4, k5, k6]β=1\displaystyle\begin{array}[]{r@{}ll}\textbf{Out:~{}~{}}&\textit{Globally}:&\texttt{[k1, k2, k3, k4, k5, k6,}\\ &&\texttt{x1(0), x2(0), x3(0),}\\ &&\texttt{x4(0), x5(0), x6(0)]}\\ &\textit{Locally not Globally}:&\texttt{[]}\\ &\textit{Non-Identifiable}:&\texttt{[]}\\ &\textit{Single-Experiment}:&\texttt{[k1, k2, k3, k4, k5, k6]}\\ &\textit{Multi-Experiment}:&\texttt{[k1, k2, k3, k4, k5, k6]}\\ &\beta=&1\end{array}