Skip to main content
  • Research article
  • Open access
  • Published:

Understanding dynamics using sensitivity analysis: caveat and solution



Parametric sensitivity analysis (PSA) has become one of the most commonly used tools in computational systems biology, in which the sensitivity coefficients are used to study the parametric dependence of biological models. As many of these models describe dynamical behaviour of biological systems, the PSA has subsequently been used to elucidate important cellular processes that regulate this dynamics. However, in this paper, we show that the PSA coefficients are not suitable in inferring the mechanisms by which dynamical behaviour arises and in fact it can even lead to incorrect conclusions.


A careful interpretation of parametric perturbations used in the PSA is presented here to explain the issue of using this analysis in inferring dynamics. In short, the PSA coefficients quantify the integrated change in the system behaviour due to persistent parametric perturbations, and thus the dynamical information of when a parameter perturbation matters is lost. To get around this issue, we present a new sensitivity analysis based on impulse perturbations on system parameters, which is named impulse parametric sensitivity analysis (iPSA). The inability of PSA and the efficacy of iPSA in revealing mechanistic information of a dynamical system are illustrated using two examples involving switch activation.


The interpretation of the PSA coefficients of dynamical systems should take into account the persistent nature of parametric perturbations involved in the derivation of this analysis. The application of PSA to identify the controlling mechanism of dynamical behaviour can be misleading. By using impulse perturbations, introduced at different times, the iPSA provides the necessary information to understand how dynamics is achieved, i.e. which parameters are essential and when they become important.


Parametric sensitivity analysis (PSA) has become a must have tool in the computational systems biologists' arsenal. In most applications of this analysis, one computes sensitivity coefficients or metrics, which generally reflect the ratios between the change in a biological model output and the perturbation on system parameters that cause this change. Depending on the magnitude of the perturbations, sensitivity analyses can be classified into local (infinitesimal perturbation) and global (finite perturbation). Regardless of these classes, the interpretation of the sensitivity metrics is intuitive; parameters with large sensitivity magnitude are deemed to be important and hence considered to be the controlling factors in the system functional regulation. Consequently, one of the common uses of PSA in systems biology is to infer the importance of cellular processes or pathways and to provide mechanistic explanations for biological behaviour [15].

On a separate note, dynamics is a prominent feature of many important biological processes (e.g., oscillations in cell cycle and circadian rhythm [6, 7], switching behaviour in programmed cell death [8], and adaptation in chemotaxis [9]). Cellular homeostatic regulation, despite the name, relies on an active dynamical response, in which orchestrated events take place in response to internal and external stimuli. Thus, understanding cellular dynamics has become a prime concern in systems biology, in which mathematical modelling coupled with quantitative analysis have been used to gain insights on the mechanisms that give rise to and control the dynamic behaviour[15]. These insights can provide the molecular targets for altering system dynamic behaviour, such as in finding treatment for diseases or in (re)engineering of cellular systems.

While there are many choices of mathematical frameworks for dynamic modelling, ordinary differential equations (ODEs) are the most commonly used modelling paradigm in systems biology and have been used to describe a wide range of biological systems. In addition, ODEs are amenable to many standard quantitative and theoretical analyses, including sensitivity analysis and bifurcation analysis, for which many off-the-shelf software packages exist that provide an integrated and user-friendly computational platform for model simulations and analyses (e.g., MATLAB [10] and XPPAUT [11]). The PSA of ODE models can be readily done using software packages such as SimBiology toolbox of MATLAB [12], PottersWheel [13], Gepasi [14], Copasi [15], JDesigner/Jarnac [16], JSim [17], BioSens [18], SBML-SAT [19], and SensSB [20]. These and other software for sensitivity analysis have been summarized in the review articles by Alves. et al.[21] and Klipp. et al.[22].

Sensitivity analysis of ODE models is well established in the science and engineering literature [2332]. In systems biology, PSA has found wide applications, such as for model calibration and identifiability, model validation and reduction, identification of bottlenecking processes, elucidation of mechanisms of complex cellular behaviour, and investigation of cellular robustness [30, 33]. A few notable examples of PSA applications in dynamic biological models include programmed cell death [3439], budding yeast cell cycle control [6], IL-6 signalling pathway [1], circadian rhythm models [7, 40, 41], and coupled MAPK and PI3K signal transduction pathway [42]. In many applications, PSA is used to generate parameter ranking based on the magnitude of sensitivity coefficients, either taken at a specific time or using consolidated sensitivity metrics, such as time-integral or average or norm of sensitivity coefficients [34, 43, 44]. The parameter ranking is subsequently used to conclude about the mechanism or property (such as robustness) of the biological system behaviour [15].

In this article, we show that the dynamical aspects of cellular functional regulation cannot be inferred from the sensitivity coefficients of PSA, neither directly nor as consolidated sensitivity metrics. More importantly, the corresponding parameter rankings from PSA can give erroneous inference about the controlling mechanisms. Briefly, the reason stems from the fact that in PSA, perturbations are introduced on system parameters, which are time-invariant or static. In other words, these parametric perturbations are persistent and their effects on the system behaviour are integrated over time. Therefore, while PSA can indicate which parameter perturbations are important, it does not point to when these perturbations matter. This problem is illustrated using local PSA of two examples: a synthetic network and a model of programmed cell death [37]. Although the illustration here was done using local sensitivity analysis, the same issue generally applies to global PSA.

To overcome this issue, a new parametric sensitivity analysis is developed in this work. This analysis differs from the classical PSA in the manner of which perturbations are introduced on model parameters, specifically using impulses, and thus is named impulse parametric sensitivity analysis (iPSA). By analyzing the consequence of impulse parameter perturbations introduced at different times, the iPSA provides time-varying, mechanistic explanation on how system dynamics is carried out. The new insights from the iPSA are demonstrated using the same two examples mentioned above.

Results and discussion

Simple network model

To illustrate the shortcoming of local PSA in analyzing system dynamics, consider a simple six state model involving three reactions with Michealis-Menten (MM) kinetics, as shown in Figure 1(a) (model parameters, rate kinetics and initial concentrations are given in Additional File 1: Supplementary Table S1). In this network, the activation of x6 followed a switch-like dynamics in response to the stimulus x1, as illustrated in Figure 1(b) (nominal). The model describes two pathways that contribute to x6 activation: a direct x2 pathway and an indirect x2, x3, and x5 pathway.

Figure 1
figure 1

A simple network model. (a) A simple network with 6 states and 4 reactions. Straight arrows connect substrate to product and the dotted arrows indicate enzymatic activity. Details of model equations and parameter values are given in Table S1 (see Additional file 1: Supplementary Table S1). (b) Activation of x 6 under a constant stimulus of x1 = 1: complete network (+), indirect pathway knock-out (*) and direct pathway knockout (x).

In this example, in silico knock-out experiments were performed by removing each pathway individually in order to assess the dominance of one pathway over the other in x6 activation. Both full network and knock-out (KO) simulations were performed under a stimulus of x1(t0 = 0) = 1. As illustrated in Figure 1(b), while the initial x6 activation in the indirect pathway KO remained the same as that of the original model, the switch-like activation was much less pronounced. On the other hand, the original switching behaviour was preserved in the direct pathway KO, but the switching time was delayed due to a slower initial activation. Taken together, these KO simulations suggested that the x6 activation is mainly accomplished through the indirect pathway, while the direct pathway contributes mainly to the initial x6 activation.

Parametric sensitivity analysis for dynamical systems: A caveat

Local parametric sensitivity analysis was also used to study the pathway dominance in this simple network. Mathematically, the parametric sensitivity coefficient is given as


where x i is the i-th state in an ODE model and p j is the j-th kinetic parameter of an ODE model (for a detailed description of sensitivity coefficient derivation, see Methods). These sensitivity coefficients describe the change in system output (state trajectory) at time t with respect to (an infinitesimal) perturbation on the system parameter values at time τ. Here, the PSA was performed for the same stimulus x1(0) = 1 with τ = 0 and the sensitivity coefficients were computed for the time range of 0-15 time units. The term local sensitivity analysis refers to the fact that the results will depend on the nominal parameter values around which the derivatives in Eq. (1) are calculated.

The sensitivities of x6 with respect to all model parameters are ranked in Figure 2 using consolidated sensitivity metrics, i.e. infinite norm [44] (Figure 2(a)), FIM [43] (Figure 2(b)), and time-integral [34] (Figure 2(c)), and using sensitivity magnitudes at switching time (t = 7.12 time units; Figure 2(d)). The conclusion from these rankings was the same: (1) the largest sensitivity was associated with the kinetics of x1 conversion to x2 and (2) the direct pathway (r2) parameters have larger sensitivities than those from the indirect pathway (r4), suggesting larger influence of the direct pathway on the x6 activation. Hence, the conclusion from the PSA is in direct contradiction with the findings from in silico KO experiments.

Figure 2
figure 2

Local parametric sensitivity analysis of x 6 activation under x 1 stimulus. The bar graphs show the consolidated sensitivity metrics of x 6 with respect to model parameters based on (a) infinite norm, (b) Fisher Information Matrix (FIM) and (c) time integrated sensitivity coefficients; and (d) the sensitivity magnitudes at switching time (t = 7.12 time units). The parameter numbers refer to the reactions as shown in Figure 1, where the subscripts f denotes the forward rate constant and k and v denotes the rate constants of Michealis-Menten kinetics.

This discrepancy can be explained by taking a closer look at the way parametric sensitivity coefficients in eqn. (1) are calculated:


where τ = t0 is the usual perturbation time, is the time derivative of sensitivity coefficient (see Methods for detail) and H(t) is the Heaviside step function. In this case, since model parameters are static or time-invariant, the parametric perturbations in the PSA consist of step changes in the parameter values, as depicted in Figure 3(a). Hence, the sensitivity coefficients at time t represent an integrated or accumulated change in the states from τ to t due to a persistent parameter change started at time τ (see Figure 3(b)). Indeed, substituting the full equation of (see Methods) in eqn. (2) gives

Figure 3
figure 3

Parametric perturbations in sensitivity analysis. An illustration of parametric perturbation and its effect on system dynamics in (a-b) PSA, (c-d) iPSA and (e-f) pulse-approximation of iPSA. Solid lines represent the nominal and the dashed lines show the perturbed trajectory, respectively. Figures are not drawn to scale.


Here, one sees two terms in the integrand that contribute to the sensitivity coefficients at time t: (1) the first is related to the (integrated) sensitivities that are carried over from the initial perturbation time τ and (2) the second accounts for the instantaneous rate changes due to the parametric perturbations that still persist at time Thus, in the PSA, a large sensitivity magnitude of Si,j (t,τ) indicates the importance of the j-th parameter in time window of τ and t, during which the perturbation is applied to the system. Hence, the use of these coefficients to infer the dynamical importance of parameters is inappropriate and can even be misleading.

For the reason above, the PSA of the simple network model gave an incorrect conclusion regarding direct versus indirect pathway activating x6. As seen in the in silico KO experiments, the direct pathway regulates the initial activation of x6, while the actual switching is carried out by the indirect pathway. In the PSA of this model (τ = 0), the early importance of the direct pathway and also the reaction r1 persisted beyond the initial times in the sensitivity coefficients due to the aforementioned integrated effect. In this case, the importance of the indirect pathway was not apparent from the parameter sensitivity rankings in the background of large (integrated) sensitivities with respect to r1 and the direct pathway. Correspondingly, the use of any time-consolidated sensitivity metrics will only worsen this problem.

Impulse parametric sensitivity analysis (iPSA)

As the problem with local PSA mentioned above is rooted from the persistent parameter perturbations, which is also done in global PSA [45], a new sensitivity analysis is formulated here that introduces impulse perturbations to model parameters as shown in Figure 3(c-d). The corresponding impulse sensitivity coefficients iSi,j (t,τ) reflect the change in the i-th state x i at time t due to an impulse perturbation on the j-th parameter p j at time τ (see Methods for the derivation and definition). Since impulse perturbations on parameters cause an immediate state changes at the perturbation time (see Methods), the inference of dynamical parametric importance can be obtained from impulse sensitivities by varying the time of perturbation. However, as with the local PSA, impulse perturbations are also local in nature and thus the impulse sensitivities will depend on the nominal parameter values. The global equivalent of iPSA can be formulated using pulse perturbations and is currently under development. Finally, from time-varying impulse perturbations, the iPSA can give answer to the following questions about state dynamics: which are the important parameters and when do they become important?

The iPSA was also applied to the simple network model under the same stimulus and for the same time range as above. Since the impulse perturbations are delivered at different times, the impulse sensitivity iSi,j (t,τ) has two-time dependence, with respect to the time of perturbation (τ) and the time of observation (t). Figure 4 shows the iPSA sensitivity coefficients of x6 at the end of simulation time (t = 15 time units), with respect to the four most important parameters at different perturbation times (for complete iPSA sensitivities, see Additional File 1: Supplementary Figure S1). In agreement with the KO simulations and in contrast to the findings from PSA, the impulse sensitivities gave support to the dominance of the indirect pathway. Specifically, the results showed that x6 activation: (1) is initiated by r1 (high initial sensitivity to parameter kf1); (2) is accomplished mainly by the direct pathway prior to the switching time (sensitivity to r 2 is higher than to r4 during these times); and (3) is subsequently carried by the indirect pathway (highest sensitivity to r4 during switching times). A higher resolution analysis using a heat map of the complete iPSA coefficients with t and τ between 0 and 15 time units gave the same conclusion (see Additional file 1: Supplementary Figure S2).

Figure 4
figure 4

Impulse parametric sensitivity analysis of x 6 activation under x 1 stimulus. iPSA coefficients of x6 with respect to impulse parametric perturbations in kf1(+), kv2(*), kv3(x) and kv4 (•). The parameter definitions and values are reported in Additional file 1: Supplementary Table S1.

Fas-induced cell death model of human Jurkat T-cells

In the second case study, we consider a more complex biological network taken from the modelling of programmed cell death in Jurkat cancer T-cells (see Figure 5). The model equations and parameters were identified from experimental data [37] (see Additional file 1: Supplementary Table S2 for detailed reaction rates and parameters). In this network, the cell death is decided by the cleaving of procaspase-3 into caspase-3 [46], in response to FasL stimulus. The model simulation showed that caspase-3 is switched ON at around t = 6000 seconds (see Figure 5, inset) and like the simple network above, there exist two activating pathways: the direct caspase-8 (type-I) and the indirect mitochondria-dependent pathway (type-II). As in the previous case study, the classical PSA and iPSA were applied to this network to elucidate the pathway dominance.

Figure 5
figure 5

Fas-induced apoptotic cell death model of human Jurkat T-cell lines [37]. Detailed model equations and parameter values are given in Additional file 1: Supplementary Table S2. The activation of cell death effector, caspase-3 follows a switch like response to a constant FasL stimulus of 2 nM (see inset).

The classical and iPSA analyses were calculated under a constant FasL stimulation (FasL = 2 nM) over the time range of 10,000 seconds. The rankings of the important parameters that control caspase-3 according to the PSA are shown in Figure 6, using consolidated metrics: infinite norm [44](Figure 6(a)), FIM [43] (Figure 6(b)), and time integral [34](Figure 6(c)), and using the sensitivity magnitudes at switching time (Figure 6(d)) (see Additional file 1: Supplementary Figure S3 for detailed sensitivity rankings). From these rankings, one could not obtain any definitive conclusion regarding the dominance of one pathway over the other. On the contrary, the iPSA sensitivity of caspase-3 in Figure 7 clearly supported a type-II dependent caspase-3 switching with an early type-I dependent activation, in agreement with two previous analyses of this model using the Green's function matrix [47] and model reduction [48]. In Figure 7, the high caspase-3 sensitivities with respect to impulse perturbations to upstream reactions were expected during the initial times, as these served as the early cell death signal response. The cleaving of procaspase-3 was carried out by caspase-8 directly (r5 of type-I) before t = 4000 seconds, after which the mitochondrial pathway (r14 of type-II) became the main route of activating caspase-3 (see Additional file 1: Supplementary Figures S4 and S5 for more detail).

Figure 6
figure 6

Local parametric sensitivity analysis of the programmed cell death model. (a-c) The bar graphs represent the ten largest sensitivity metrics of caspase-3 based on infinite norm, Fisher Information Matrix (FIM), time integrated sensitivity coefficients, and sensitivity magnitudes at the switching time (t = 6060s). The complete parameter ranking is given in Additional file 1: Supplementary Figure S3. The parameter numbers refer to the reactions shown in Figure 5, where the subscripts f and r denote forward and backward rate constants for reversible reactions and k denotes the rate constants for irreversible reactions.

Figure 7
figure 7

Impulse parametric sensitivity analysis of the programmed cell death model. The iPSA sensitivity of caspase-3 with respect to impulse parametric perturbations in J1_f(+), J5_k(*), J9_f(x) and J14_k(•). The complete iPSA sensitivities are given in Additional file 1: Supplementary Figures S4 and S5. The parameter numbers refer to the reactions shown in Figure 6, where the subscripts f and r denote forward and backward rate constants for reversible reactions and k denotes the rate constants for irreversible reactions.

The discrepancy between the PSA and iPSA results can again be explained in the context of persistent versus impulse perturbations. As seen in the PSA parameter rankings in Figure 6 and following the insights offered by the iPSA in Figure 7, the effect of perturbing early response processes, including type-I pathway, was integrated over time in the PSA. For example, the highest ranked parameters in the PSA were associated with the first three reactions, r1 to r3. Such integration masked the dynamical importance of different parameters in this analysis. The conclusion from the iPSA is in agreement with the simulations of KOs of type-I and type-II pathways in Figure 8. While type-II knock-out could describe the caspase-3 activation at early times, it failed to capture the switching behaviour. On the other hand, the type-I knock-out was able to reproduce the switching of caspase-3, albeit with a short delay.

Figure 8
figure 8

In silico knock-out of the programmed cell death network. The activation of caspase-3 are simulated for nominal network (-) and reduced networks by knocking out type-I (•) and type-II pathway (x).

The two examples above illustrate the problem of using the classical PSA in identifying the controlling mechanisms of a dynamical system. Of course, this does not mean that the PSA of dynamical models is incorrect, but rather the interpretation of the sensitivity coefficients should be carefully managed. In particular, a large sensitivity magnitude with respect to a parameter suggests the importance of this parameter in the time period between the perturbation time τ and the state observation time t. In contrast, the iPSA is developed with dynamics in mind, where the impact of a single perturbation on the system is realized only at the perturbation time and subsequently they are delivered at varying perturbation times. By doing so, the iPSA coefficients can elucidate the way system dynamics x(t) is achieved, by indicating which and when parameters or processes are essential. Because of the persistent nature of perturbations used in the PSA, it is still not possible to reproduce the conclusions of the iPSA by varying the time of perturbation (see Additional file 1: Supplementary Figures S6 and S7).


While classical parametric sensitivity analysis provides a powerful tool to understand the parametric dependence of biological behaviour, its suitability in inferring mechanisms of dynamic behaviour has not been properly addressed. The two case studies here illustrated the caveat of using local PSA for such purpose. The issue mainly arose from the information needed to do this inference, where one needs to know not only which parameters are critical, but also when they matter. However, the persistent parametric perturbations in standard PSA are incapable of providing this information as the sensitivity coefficients represent an integrated effect. A new sensitivity analysis, called impulse parametric sensitivity analysis (iPSA), was developed with dynamical systems in mind. In particular, the iPSA makes use of local impulse perturbations introduced at different times to produce the necessary information for understanding dynamics. The application of iPSA to the case studies was able to correctly pinpoint the mechanisms responsible for dynamical system behaviour, while local PSA failed in these cases. Since the discrepancy between PSA and iPSA arises from a fundamental difference in the manner of which parametric perturbations are realized (i.e. persistent vs. impulse), the same caveat and solution can be generalized to the global PSA, in which the perturbations are no longer infinitesimal.


Mathematical Models

Ordinary differential equation models of dynamic systems can generally be written as:


In biological models, the state x nis typically the concentration vector of biomolecular species, such as mRNAs and proteins, while the function f is the constitutive, often nonlinear, rate equation. The right hand side of the ODE captures the generation and consumption of biomolecules due to a variety of processes in the cell (e.g. transcription, translation, phosphorylation and dephosphorylation, etc), the rates of which depend on a set of kinetic parameters that are consolidated in the vector . Since the initial conditions x0 can be treated in the same way as model parameters, the aggregate vector p m+nis used here to denote the combined parameters and initial conditions, i.e. .

Local parametric sensitivity analysis (PSA)

The effect of parameter perturbations can be written in a Taylor series expansion:


where the partial derivatives ∂x i /∂p j 's are the first-order sensitivity coefficients, describing the linear change in the state x i at any time t with respect to an infinitesimal perturbation to the parameter p j . In general, the parametric perturbation can be introduced at any time τ(t0τt) [27], i.e.


but in the PSA, the perturbation time τ is commonly taken to be the initial time t0. Hence, the argument τ is typically dropped out of eqn. (6) and the sensitivity coefficients only carry a single time dependence on the observation time t[13, 5, 30, 44, 49, 50]. The higher order sensitivity coefficients in the Taylor series expansion are less commonly computed, and hence the focus of the current work is only on the first-order sensitivities. Because the magnitude of perturbations are infinitesimally small, the sensitivity coefficients will depend on the nominal or baseline parameter values, and thus the classical PSA is considered a local analysis.

The sensitivity coefficients of an ODE model can be computed by directly differentiating the model eqn. (4) with respect to parameters p, giving the differential equation for the sensitivity coefficients as:


where ∂f/∂x is also known as the Jacobian matrix and 0n×mand In×nare n × m zero and n × n identity matrices, respectively. The parametric sensitivity coefficients in eqn. (7) need to be solved simultaneously with the ODE model in eqn. (4), which is called the direct method [28]. As the state and parameter values may span a large range of magnitudes, normalized sensitivity values are often used to compare among states and parameters and to generate parameter ranking, which is given by:


In most systems biology applications of PSA, parameter rankings are generated from the sensitivity coefficients, either directly or using some consolidated sensitivity metrics. In this article, four sensitivity metrics are used to rank parameters based on the PSA results:


where the indices i and j again denote the i-th state and j-th parameter, and Sinf, SFIM, Sint and are the sensitivity metrics based on infinite norm [44], Fisher information matrix [43], time integral [34] and sensitivity magnitude at a particular time, respectively.

Impulse parametric sensitivity analysis (iPSA)

The derivation of the iPSA coefficient follows the illustration in Figure 3(e). The sensitivity coefficients of iPSA are constructed by quantifying the ratio between the change in the state x i at time t and the causative pulse perturbation of size Δp j τ for a duration of Δt, which is applied to the parameter p j at time τ, in the limit Δp j and Δτ tending to zero. The first step of the derivation is to quantify the change in all states x at the end of the pulse perturbation, i.e. at time τ + Δτ, using the Taylor series expansion:


In the next step, the change Δx(τ + Δτ) is translated to the change in the state x i at time t using the Green's function matrix (GFM) Sx(t,τ + Δτ) [47]. The (i,j)-th element of the GFM represents the sensitivity of the state x i at time t to alteration in the state x j at some previous time τ, i.e.


Thus, the change Δx i (t) due to the pulse perturbation is given by


Subsequently, substitution of eqn. (10) in eqn. (12) gives


Then, taking Taylor series expansion of the parametric sensitivities around the time τ and dividing both sides by Δp j , one arrives with:


Finally, taking the limit as Δp j τ → 0 such that the pulse perturbation becomes an impulse, the iPSA coefficient is obtained as:


Note that at t = τ, the iPSA coefficient reduces to


since Sx(τ,τ) = I. In other words, the impulse parameter perturbation causes an immediate change in the state x at time τ. By rewriting eqn. (15) as:


one can further see that the impact of this impulse perturbation takes effect only at the perturbation time τ and that the consequence on the state trajectory is equivalent to perturbing the states themselves, similar to the GFM analysis. Like in the PSA, the iPSA coefficients should be normalized for comparison and parameter ranking purposes, according to:






Ordinary Differential Equations


Parametric Sensitivity Analysis


Impulse Parametric Sensitivity Analysis.


  1. Chu Y, Jayaraman A, Hahn J: Parameter sensitivity analysis of IL-6 signalling pathways. IET Syst Biol 2007, 1: 342-352. 10.1049/iet-syb:20060053

    Article  CAS  PubMed  Google Scholar 

  2. Ihekwaba AE, Broomhead DS, Grimley RL, Benson N, Kell DB: Sensitivity analysis of parameters controlling oscillatory signalling in the NF-kappaB pathway: the roles of IKK and IkappaBalpha. Syst Biol (Stevenage) 2004, 1: 93-103. 10.1049/sb:20045009

    Article  CAS  Google Scholar 

  3. Ihekwaba AE, Broomhead DS, Grimley R, Benson N, White MR, Kell DB: Synergistic control of oscillations in the NF-kappaB signalling pathway. Syst Biol (Stevenage) 2005, 152: 153-160.

    Article  CAS  Google Scholar 

  4. Adler P, Peterson H, Agius P, Reimand J, Vilo J: Ranking Genes by Their Co-expression to Subsets of Pathway Members. Annals of the New York Academy of Sciences 2009, 1158: 1-13. 10.1111/j.1749-6632.2008.03747.x

    Article  CAS  PubMed  Google Scholar 

  5. Iwamoto K, Tashima Y, Hamada H, Eguchi Y, Okamoto M: Mathematical modeling and sensitivity analysis of G1/S phase in the cell cycle including the DNA-damage signal transduction pathway. Biosystems 2008, 94: 109-117. 10.1016/j.biosystems.2008.05.016

    Article  CAS  PubMed  Google Scholar 

  6. Lovrics A, Zsély GyI, Csikász-Nagy A, Zádor J, Turányi T, Novák B: Analysis of a budding yeast cell cycle model using the shapes of local sensitivity functions. International Journal of Chemical Kinetics 2008, 40: 710-720. 10.1002/kin.20366

    Article  CAS  Google Scholar 

  7. Jin Y, Peng X, Liang Y, Ma J: Uniform design-based sensitivity analysis of circadian rhythm model in Neurospora. Computers and Chemical Engineering 2008, 32: 1956-1962. 10.1016/j.compchemeng.2007.10.013

    Article  CAS  Google Scholar 

  8. Zhang T, Brazhnik P, Tyson JJ: Computational analysis of dynamical responses to the intrinsic pathway of programmed cell death. Biophys J 2009, 97: 415-434. 10.1016/j.bpj.2009.04.053

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  9. Hansen CH, Endres RG, Wingreen NS: Chemotaxis in Escherichia coli: a molecular model for robust precise adaptation. PLoS Comput Biol 2008, 4: e1. 10.1371/journal.pcbi.0040001

    Article  PubMed Central  PubMed  Google Scholar 

  10. MATLAB

  11. XPPAUT

  12. Schmidt H, Jirstrand M: Systems Biology Toolbox for MATLAB: a computational platform for research in systems biology. Bioinformatics 2006, 22: 514-515. 10.1093/bioinformatics/bti799

    Article  CAS  PubMed  Google Scholar 

  13. Maiwald T, Timmer J: Dynamical modeling and multi-experiment fitting with PottersWheel. Bioinformatics 2008, 24: 2037-2043. 10.1093/bioinformatics/btn350

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  14. Mendes P: GEPASI: a software package for modelling the dynamics, steady states and control of biochemical and other systems. Comput Appl Biosci 1993, 9: 563-571.

    CAS  PubMed  Google Scholar 

  15. Hoops S, Sahle S, Gauges R, Lee C, Pahle J, Simus N, Singhal M, Xu L, Mendes P, Kummer U: COPASI--a COmplex PAthway SImulator. Bioinformatics 2006, 22: 3067-3074. 10.1093/bioinformatics/btl485

    Article  CAS  PubMed  Google Scholar 

  16. JDesigner/Jarnac

  17. JSim

  18. BioSens

  19. Zi Z, Zheng Y, Rundell AE, Klipp E: SBML-SAT: a systems biology markup language (SBML) based sensitivity analysis tool. BMC Bioinformatics 2008, 9: 342. 10.1186/1471-2105-9-342

    Article  PubMed Central  PubMed  Google Scholar 

  20. Rodriguez-Fernandez M, Banga JR: SensSB: a software toolbox for the development and sensitivity analysis of systems biology models. Bioinformatics 26: 1675-1676. 10.1093/bioinformatics/btq242

  21. Alves R, Antunes F, Salvador A: Tools for kinetic modeling of biochemical networks. Nat Biotechnol 2006, 24: 667-672. 10.1038/nbt0606-667

    Article  CAS  PubMed  Google Scholar 

  22. Klipp E, Liebermeister W, Helbig A, Kowald A, Schaber J: Systems biology standards--the community speaks. Nat Biotechnol 2007, 25: 390-391. 10.1038/nbt0407-390

    Article  CAS  PubMed  Google Scholar 

  23. Saltelli A, Chan K, Scott EM: Sensitivity Analysis: Gauging the Worth of Scientific Models. John Wiley \& Sons, Ltd.; 2000.

    Google Scholar 

  24. Saltelli A, Ratto M, Tarantola S, Campolongo F: Sensitivity analysis for chemical models. Chem Rev 2005, 105: 2811-2827. 10.1021/cr040659d

    Article  CAS  PubMed  Google Scholar 

  25. Saltelli A, Tarantola S, Campolongo F: Sensitivity Analysis as an Ingredient of Modeling. Stat Sci 2000, 15: 377-395. 10.1214/ss/1009213004

    Article  Google Scholar 

  26. Saltelli A, Tarantola S, Campolongo F, Ratto M: Sensitivity Analysis in Practice: A Guide to Assessing Scientific Models. John Wiley \& Sons, Ltd.; 2004.

    Google Scholar 

  27. Turányi T: Sensitivity analysis of complex kinetic systems. Tools and applications. Journal of Mathematical Chemistry 1990, 5: 203-248.

    Article  Google Scholar 

  28. Varma A, Morbidelli M, Wu H: Parametirc Sensitivity in Chemical Systems. Cambridge University Press, Cambridge, UK; 1999.

    Book  Google Scholar 

  29. Frey HC, Patil SR: Identification and review of sensitivity analysis methods. Risk Anal 2002, 22: 553-578. 10.1111/0272-4332.00039

    Article  PubMed  Google Scholar 

  30. Ingalls B: Sensitivity analysis: from model parameters to system behaviour. Essays Biochem 2008, 45: 177-193. 10.1042/BSE0450177

    Article  CAS  PubMed  Google Scholar 

  31. Marino S, Hogue IB, Ray CJ, Kirschner DE: A methodology for performing global uncertainty and sensitivity analysis in systems biology. J Theor Biol 2008, 254: 178-196. 10.1016/j.jtbi.2008.04.011

    Article  PubMed Central  PubMed  Google Scholar 

  32. Rabitz H, Kramer M, Dacol D: Sensitivity analysis in chemical kinetics. Annual review of Physical Chemistry 1983, 34: 419-461. 10.1146/annurev.pc.34.100183.002223

    Article  CAS  Google Scholar 

  33. Stelling Jr, Sauer U, Szallasi Z, Doyle FJ, Doyle J: Robustness of cellular functions. Cell 2004, 118: 675-685. 10.1016/j.cell.2004.09.008

    Article  CAS  PubMed  Google Scholar 

  34. Bentele M, Lavrik I, Ulrich M, Stober S, Heermann DW, Kalthoff H, Krammer PH, Eils R: Mathematical modeling reveals threshold mechanism in CD95-induced apoptosis. J Cell Biol 2004, 166: 839-851. 10.1083/jcb.200404158

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  35. Eissing T, Allgower F, Bullinger E: Robustness properties of apoptosis models with respect to parameter variations and intrinsic noise. Syst Biol (Stevenage) 2005, 152: 221-228.

    Article  CAS  Google Scholar 

  36. Hua F, Cornejo MG, Cardone MH, Stokes CL, Lauffenburger DA: Effects of Bcl-2 levels on Fas signaling-induced caspase-3 activation: molecular genetic tests of computational model predictions. J Immunol 2005, 175: 985-995.

    Article  CAS  PubMed  Google Scholar 

  37. Hua F, Hautaniemi S, Yokoo R, Lauffenburger DA: Integrated mechanistic and data-driven modelling for multivariate analysis of signalling pathways. J R Soc Interface 2006, 3: 515-526. 10.1098/rsif.2005.0109

    Article  PubMed Central  PubMed  Google Scholar 

  38. Shoemaker JE, Doyle III FJ: Identifying Fragilities in Biochemical Networks: Robust Performance Analysis of Fas Signaling-Induced Apoptosis. Biophys J 2008.

    Google Scholar 

  39. Aldridge BB, Haller G, Sorger PK, Lauffenburger DA: Direct Lyapunov exponent analysis enables parametric study of transient signalling governing cell behaviour. Syst Biol (Stevenage) 2006, 153: 425-432.

    Article  CAS  Google Scholar 

  40. Stelling Jr, Gilles ED, Doyle FJ: Robustness properties of circadian clock architectures. Proc Natl Acad Sci USA 2004, 101: 13210-13215. 10.1073/pnas.0401463101

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  41. Gunawan R, Doyle FJ: Phase sensitivity analysis of circadian rhythm entrainment. J Biol Rhythms 2007, 22: 180-194. 10.1177/0748730407299194

    Article  PubMed  Google Scholar 

  42. Hu D, Yuan JM: Time-dependent sensitivity analysis of biological networks: coupled MAPK and PI3K signal transduction pathways. J Phys Chem A 2006, 110: 5361-5370. 10.1021/jp0561975

    Article  CAS  PubMed  Google Scholar 

  43. Yue H, Brown M, He F, Jia J, Kell DB: Sensitivity analysis and robust experimental design of a signal transduction pathway system. International Journal of Chemical Kinetics 2008, 40: 730-741. 10.1002/kin.20369

    Article  CAS  Google Scholar 

  44. Yue H, Brown M, Knowles J, Wang H, Broomhead DS, Kell DB: Insights into the behaviour of systems biology models from dynamic sensitivity and identifiability analysis: a case study of an NF-kappaB signalling pathway. Mol Biosyst 2006, 2: 640-649. 10.1039/b609442b

    Article  CAS  PubMed  Google Scholar 

  45. Saltelli A: Global sensitivity analysis: the primer. Chichester, England: John Wiley; 2008.

    Google Scholar 

  46. Reed JC, Doctor KS, Godzik A: The domains of apoptosis: a genomics perspective. Sci STKE 2004 2004, re9. 10.1126/stke.2392004re9

    Google Scholar 

  47. Perumal TM, Wu Y, Gunawan R: Dynamical analysis of cellular networks based on the Green's function matrix. J Theor Biol 2009, 261: 248-259. 10.1016/j.jtbi.2009.07.037

    Article  PubMed  Google Scholar 

  48. Okazaki N, Asano R, Kinoshita T, Chuman H: Simple computational models of type I/type II cells in Fas signaling-induced apoptosis. J Theor Biol 2008, 250: 621-633. 10.1016/j.jtbi.2007.10.030

    Article  CAS  PubMed  Google Scholar 

  49. Yamada S, Shiono S, Joo A, Yoshimura A: Control mechanism of JAK/STAT signal transduction pathway. FEBS Lett 2003, 534: 190-196. 10.1016/S0014-5793(02)03842-5

    Article  CAS  PubMed  Google Scholar 

  50. Zi Z, Cho KH, Sung MH, Xia X, Zheng J, Sun Z: In silico identification of the key components and steps in IFN-gamma induced JAK-STAT signaling pathway. FEBS Lett 2005, 579: 1101-1108. 10.1016/j.febslet.2005.01.009

    Article  CAS  PubMed  Google Scholar 

Download references


TMP is supported by the Singapore Millennium Foundation scholarship.

Author information

Authors and Affiliations


Corresponding author

Correspondence to Rudiyanto Gunawan.

Additional information

Authors' contributions

TMP performed the experiments and wrote the manuscript. RG supervised the study and wrote the manuscript. All authors have read and approved the final manuscript.

Electronic supplementary material


Additional file 1:Supplementary Material. Detailed results of iPSA and PSA, including model equations and parameters, of the simple network model and Fas-induced cell death model (PDF 574 KB)

Authors’ original submitted files for images

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and permissions

About this article

Cite this article

Perumal, T.M., Gunawan, R. Understanding dynamics using sensitivity analysis: caveat and solution. BMC Syst Biol 5, 41 (2011).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: