Open access peer-reviewed chapter

Sensitivity Analysis: A Useful Tool for Bifurcation Analysis

By Raheem Gul and Stefan Bernhard

Submitted: April 19th 2017Reviewed: November 9th 2017Published: May 2nd 2018

DOI: 10.5772/intechopen.72345

Abstract

Sensitivity analysis and bifurcation analysis are closely related to each other. In sensitivity analysis, especially global sensitivity analysis the effects of input parameter spaces on output quantities of interest are studied. On the other hand, in bifurcation analysis the critical points within feasible regions of parameters are detected where the long-term dynamics changes qualitatively. Prior to bifurcation analysis, it is important to identify the bifurcation parameters. In complex and computationally expensive problems which consist plenty of uncertain parameters, it is essential to find a set of bifurcation parameters before bifurcation analysis. Global sensitivity analysis is a powerful tool to identify the bifurcation parameters which contribute most on output uncertainty. Global sensitivity analysis is the first step toward bifurcation analysis which helps in dimension reduction during bifurcation analysis. As an example, in this chapter, a multi compartment, lumped-parameter model of an arm artery is considered and global sensitivity analysis (Sobol’s method) is applied to identify the bifurcation parameters of the arm arteries.

Keywords

• lumped parameter model
• arm arteries
• sensitivity analysis
• bifurcation analysis
• bifurcation parameters
• Sobol’s method

1. Introduction

Sensitivity analysis and bifurcation analysis are closely related to each other. In sensitivity analysis, we study how the uncertainty in the output of a mathematical model or system (numerical or otherwise) can be apportioned to different sources of uncertainty in its inputs [1]. On the other hand, in bifurcation analysis the critical points within the feasible regions of parameters are detected where the long-term dynamics changes qualitatively [2]. Prior to the bifurcation analysis, it is important to identify the bifurcation parameters in complex and computationally expensive problems that consist plenty of uncertain parameters. Sensitivity analysis is a powerful tool to identify the bifurcation parameters which contribute most on output uncertainty. Also, sensitivity analysis helps in dimension reduction during the bifurcation analysis by fixing less influential parameters on their nominal values.

Sensitivity analysis can be divided into two categories, local sensitivity analysis (LSA) and global sensitivity analysis (GSA). In LSA a parameter value is perturbed around its nominal values at a time, keeping other parameters fixed on their nominal values. The procedure is repeated for all parameters one by one to study their impact on output variables. LSA techniques are simple, easy to implement and computationally less expensive. On the other hand, LSA is not suitable for non-linear models and does not explore the impact of entire parameter spaces and their interactions effects on output variables [3, 4]. In order to overcome the limitations associated with the LSA, GSA can be used. In GSA, the analysis is performed over entire feasible regions of the input parameters and quantifies the impact of parameter interactions on output variables. The only deficiency related to the GSA is its computational cost [5, 6, 7, 8, 9, 10, 11, 12] (Figure 1).

In this chapter, the main questions of interest are:

1. How to identify the bifurcation parameters in a model having plenty of input parameters?

2. Which parameters could be exempted from the bifurcation analysis (dimension reduction)?

This chapter seeks to answers the above-mentioned questions using a simplified 5-steps procedure of uncertainty and sensitivity analysis. As an example, a multi-compartment, lumped-parameter model of arm arteries is considered [4] and global sensitivity analysis (Sobol’s method) is applied to identify the bifurcation parameters (electrical) of the arm arteries.

2. Lumped-parameter model of the arm arteries

In this section, the major arteries of the arm are divided in to number of non-terminal and terminal arterial segments (nodes). The total number of arterial segments, Ns=15including 12 non-terminal and 3 terminal segments. Each non-terminal and terminal arterial segment is represented by its corresponding non-terminal and terminal electrical circuit.

Applying Kirchhoff’s current and voltage laws on electrical representation of arm arteries, the following mathematical equations for pressure and flow are obtained:

Pressure and flow equations at non-terminal nodes:

Flow equation:

q˙i=pi1piRiqiLi,i=1,2,3,15andi11,13q˙11=p6p11R11q11L11q˙13=p11p13R13q13L13E1

Pressure equation:

p˙i=qiqi+1Ci,i=1,2,3,15andi6,11p˙6=q6q11q7C6,p˙11=q11q12q13C11atbifurcationE2

Pressure and flow equations at terminal nodes:

q˙in=2pin2piRiqinLip˙i=qinqoutCiq˙out=2pi2pout2RbqoutLi,i=10,12,15E3

where, Ri,CiandLiis the blood flow resistance, compliance of the vessel and blood inertia of ithsegment of the arm arteries respectively. The electrical parameters RiCiLiof ithsegments are related with structural parameters Eilidihias,

Ri=8νliπd24,Ci=ρliπd22,Li=2πd22liEiliE4

where, Eiis the Young modulus, lidenotes length of the vessel, dis the diameter of the vessel and hirepresents the wall thickness of ithsegment of the vessel. Moreover, ν(0.004 Pa s) is the blood viscosity and ρ(1050 kgm3) is the blood density. The nominal values of all parameters of arm segments are given in Table 1. The geometry along with the values of the parameters is taken from [13, 14].

NodesEldhRCL
unitskgm2s2×105m×102m×103m×104kgs1m4×106kg1s2m4×1011kgm4×106
146.17.286.23.5397.4541.539
245.66.285.75.8684.7781.898
346.35.645.510.154.0352.648
446.35.325.312.823.5142.976
546.355.216.432.9743.369
644.64.72515.101.92.76
787.13.484.478.900.6677.838
887.13.244.31050.5319.042
987.134.2142.90.44810.55
10822.844.155.110.12073.647
11824.34.931.941.0674.844
12166.71.822.811730.083431.88
1387.94.064.740.190.93665.434
1486.73.484.650.220.806.075
1583.73.664.533.600.39583.693

Table 1.

Numerical values of parameters for each node of the arm arteries (shown in Figure 2).

The value of boundary resistance Rbon three terminal nodes is 3.24×109 kgs1m4, ν=0.004 kgs1m1and ρ=1050 kgm3[4, 5, 13, 14].

3. Uncertainty and sensitivity analysis

Uncertainty analysis (UA) and sensitivity analysis (SA) are closely related; however they represent two different disciplines. Uncertainty analysis assesses the uncertainty in model outputs caused by uncertainty of its inputs. Whereas, sensitivity analysis study the impact of input quantities of interest (QoI) on output quantities of interest (QoI). In this study, the input (QoI) are electrical parameters (Ri,Ci,Li) and output (QoI) are pressure and flow at each node of the arm arteries. Further, for uncertainty analysis Latin hypercube sampling (LHS) is used and variance-decomposition method (Sobol’s method) is used for global sensitivity analysis (GSA).

Compared to the high-dimensional cardiovascular models (3D, 2D, 1D), lumped-parameter models of the cardiovascular system (CVS) are computationally less expensive, therefore they are suitable for GSA. In our previous studies, we found that for lumped-parameter models of the CVS, the Sobol’s method is computationally less expensive as compared to the other variance-decomposition methods, like FAST and sparse grid stochastic collocation method based on Smolyak algorithm [5, 12].

3.1. The method of Sobol

The method of Sobol is the variance-decomposition method used for global sensitivity analysis. The method decomposes the output variance of a system or model into fractions and assigns them to the inputs factors. For example, given a model of the form Y=fX=fx1x2xk, where Xis the vector of Kuncertain parameters, which are independently generated within a unit hypercube i.e. xi01kfor i=1,2,3,,K. Compared to the other GSA methods, the Sobol’s method is one of the most commonly used variance-decomposition method, because of its ease of implementation. The method is primarily based on the decomposition of output Yinto summands of elementary functions in terms of increasing dimensionality [1, 8],

fx1x2xk=f0+ikfixi+iki<jkfijxixj++f1,2,3,,kx1x2x3xkE5

In Eq. (5), fis integrable, f0is a constant, fiis a function of xi, fijis a function of xiand xjand so on. Furthermore, all the terms in the functional decomposition are orthogonal, which leads toward the following definitions of the terms of the functional decomposition in term of conditional expected values.

f0=EYfixi=ExiYxif0fijxixj=ExijYxixjf0fifjE6

where, Edescribes the mathematical expectation and xidenotes all parameters except xiand so on. The total unconditional variance can be obtained by,

V=ΩKf2Xdxf02E7

From Eq. (7), the total unconditional variance can be decomposed in a similar manner like in Eq. (5) as,

V=ikVixi+iki<jkVijxixj++f1,2,3,,kx1x2x3xkE8

where, Vis the variance operator. The relationship between functions and partial variance are given by,

Vi=VxiExiYxi=VfixiVij=Vxi,xjExijYxixjViV=VfijxixjE9

Dividing both sides of the Eq. (8) by V, we get:

1=ikSixi+iki<jkSijxixj++S1,2,3,,Kx1x2x3xKE10

Where,

Si=VijV, andSij=VijVE11

where, Siis the main effect (first order sensitivity index) of the ithparameter on output uncertainty and Sijis the interaction effect of ithand jthparameters on output uncertainty. Further, the total sensitivity index, STican be calculated as,

STi=ExiVxiYxiV=1VxiExiYxiVE12

In general, the main effect is used identify the most influential parameters (bifurcation parameters) and the total effect is taken into account for those parameters which are exempted from bifurcation analysis (factor fixing). The total effect, STiof the ithparameter means main effect plus higher-order effect due to interactions of the ithparameter. In this study, the interaction effects of parameters on the output (QoI) are negligible, therefore the main effects are used for factor fixing and ranking of bifurcation parameters.

3.2. Algorithm to compute sensitivity indices

In this section, a detailed working algorithm is presented to compute the main effect, Siusing the Monte Carlo simulations, we follow the steps, given in [1, 15].

1. Generate a random numbers matrix of row dimension 2Kand column length N(the sample size) and split into two independent sampling matrices, ANKand BNKby using LHS. Where, Kis the number of uncertain model parameters.

ANK=x11x12x1Kx21x22x2KxN1xN2xNKE13
BNK=x1K+1x1K+2x12Kx2K+1x2K+2x22KxNK+1xNK+2xN2KE14

• Define matrix Ci, which is matrix A except the ithcolumn of matrix B.

CiNK=x11x12x1K+ix1Kx21x22x2K+ix2KxN1xN2xNK+ixNKE15

• Compute and save model runs for all parameter spaces using matrices A, Band Cii.e. YAtTsN=fA,YBtTsN=fBandYCitTsNK=fCi, where, tare the time points for one heart beat with period tp=0.8s, Tsrepresents the state variables (pressure and flow time series at six locations of arm artery NTs=15and Nis the total number of model runs N=4000.

• For the time dependent model outputs, we compute the time dependent main sensitivity index, of each parameter at each time-point of the pressure and flow waves, using the estimator offered by Jansen [15, 16, 17].

Sti=ViV=VxiExiYxiV=V12Nn=1NYBnYCin2V=112N×Vn=1NYBnYCin2E16

where,

V=1Nn=1NYBn2E2E17

and

E=1Nn=1NYBn2E18

The total variance Vand the expectation Eare also calculated at each time-point of pressure and flow waves with respect to each parameter.

• Finally, the main effect, Si of each parameter on the state variables is calculated.

Si=1NTs1Ntj=1NTst=0NtStitjt,i=1,2,,KE19

In Eq. (19), NTsis the number of output variables (pressure and flow time series at all locations) and Ntis the number of time-points [12].

• 3.3. Input parameters distribution

The results of the UA and SA are greatly affected by the choice of input parameters distributions. In principle, the parameters distributions should be estimated using medical data. Unfortunately, the medical data is not easy to obtained. The input parameters distributions could be chosen according to the expert opinion or using the data from the literature. Due to limited data availability, here in this work the input parameters are randomized within ±10% range of their base (nominal) values using Latin hypercube sampling (LHS).

3.4. Convergence of sensitivity indices

The method of Sobol requires NK+2number of model simulations to compute Si. The main effect, Siis computed for N= [500, 1000, 2000, 3000, 4000] model runs. It is observed that, when the total number of simulations run Nincreases from 3000 then the sensitivity indices (Si) become stable [18]. Therefore, the minimum number of simulations for each parameter to achieve convergence of sensitivity indices is around 3000.

4. Results and discussion

In this section, the sensitivity results based on main effect Siare presented. In order to calculate sensitivity time series, the method of Sobol is applied on each time point of the output QoI i.e. pressure and flow waves at each location of the arm arteries. For each parameter, there are two sensitivity time series at each segment of the arm arteries, one for the pressure and one for the flow. In total, K×NTs=45×33=1485sensitivity time series are obtained. In order to represent the sensitivity results in a compact way, mean absolute values of each pressure and flow sensitivity time series per parameter is taken. In this way, a matrix of dimension 45×33is acquired, where each entry of the matrix represents the mean absolute values of pressure and flow sensitivity time series per parameter, see Figure 3. The numbers in the boxes show the impact (%) on the output (pressure and flow) when input parameters RiCiLiare randomized within the feasible ranges of ±10%. The parameters having main effect, Si>10%on output QoI are not shown in the Figure 3. Each row in Figure 3 represents the ranking of influential (bifurcation) parameters. For convenience, the electrical parameters RiCiLi,i=1,2,3,,15that have impact greater than 10% on pressure and flow are considered as bifurcation parameters which further can be used in bifurcation analysis. For example, for pressure at node-2, L1and L2are the bifurcation parameters, see in Figure 4 (top). Whereas, for flow at node -2, L1and L2are considered as bifurcation parameters, see Figure 4 (bottom).

In a similar fashion, each row of Figure 3 represents the ranking of bifurcation parameters which further can be used in bifurcation analysis. The parameters which have main effect Si<10%can be exempted from the bifurcation analysis. The criteria for factor fixing vary from problem to problem.

5. Conclusion

In this chapter, a 5-step procedure of global sensitivity analysis is presented to identify the bifurcation parameters in a lumped-parameter model of the arm arteries. Moreover, the proposed procedure can be applied on any morphology or structure of the systemic circulation (carotid bifurcation, aorta or complete systemic circulation). The results of sensitivity analysis are useful to identify and rank the bifurcation parameters, as well as help which parameters could be exempted from the bifurcation analysis. In this particular example of the arm arteries, 23 out of 45 parameters can be excluded from the bifurcation analysis. Whereas, 22 identified as bifurcation parameters, which further can be used/studied in the bifurcation analysis.

More

© 2018 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution 3.0 License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

How to cite and reference

Cite this chapter Copy to clipboard

Raheem Gul and Stefan Bernhard (May 2nd 2018). Sensitivity Analysis: A Useful Tool for Bifurcation Analysis, Complexity in Biological and Physical Systems - Bifurcations, Solitons and Fractals, Ricardo López-Ruiz, IntechOpen, DOI: 10.5772/intechopen.72345. Available from:

Next chapter