Robert Rapadamnaba
email: robert.rapadamnaba@umontpellier.fr

Mathieu Ribatet

Bijan Mohammadi

Global sensitivity analysis for assessing the parameters importance and setting a stopping criterion in a biomedical inverse problem

Keywords: sensitivity analysis, Sobol' sensitivity indices, uncertainty quantification, hemodynamic inverse problems, parameter estimation, convergence stopping criterion

This paper shows how to obtain in addition to the standard deviations available after a data assimilation procedure based on the ensemble Kalman filter, an apportioning of the total uncertainty in the outputs of a patient-specific blood flow model into small portions of uncertainty due to input parameters. Statistical indicators generally used for identifying the importance of numerical parameters, namely the Sobol' first order and total indices, are introduced and discussed. These allow the identification of the importance rank of the different input parameters for the patient-specific blood flow model, as well as the influence of the interactions between these parameters on the model output variance. The results show that knowing the importance rank of the model input parameters during the assimilation procedure is useful to avoid unnecessary over-solving and to find a suitable stopping criterion in clinical situations where faster diagnosis is always requested. Indeed, the work permits to reduce typically by a factor of six the time to solution and most importantly with very limited extra calculation using already available information.

INTRODUCTION

Models have always been invaluable tools used to explain a phenomenon or a system, to study the effects of their different components, to communicate knowledge, and to make decision or predictions about the future behavior of a real-world phenomenon [START_REF] Schichl | Models and the History of Modeling[END_REF]. Of all models, mathematical models, which describe a given system by a set of variables and a set of equations establishing relationships between the variables, are of great importance in the research world. Indeed, most often providing a simplified version of something that is real, mathematical models are inseparable from most of scientific research disciplines to some extent [START_REF] Crosetto | Sensitivity and uncertainty analysis in spatial modelling based on gis[END_REF]3]. For instance, in many scientific fields, mathematical modeling is the starting point for analyzing a system to be understood, controlled or optimized.

Unfortunately, a mathematical model can have or include significant spurious artifacts, such as random events, measurement errors, more-or-less uncertain parameters, missing information, unrealistic nonlinearities, and also conceptual uncertainty, i.e., poor knowledge or partial, erroneous understanding of the driving forces and mechanisms in the model, or uncertainty in the model structures, assumptions and specifications [START_REF] Crosetto | Sensitivity and uncertainty analysis in spatial modelling based on gis[END_REF][START_REF] Pannell | Sensitivity analysis of normative economic models: theoretical framework and practical strategies[END_REF][START_REF] Hill | Practical use of computationally frugal model analysis methods[END_REF][START_REF] Kiureghian | Aleatory or epistemic? does it matter?[END_REF]. Moreover, it can particularly be highly complex, and as a result, the relationships between its inputs and its outputs may be poorly understood. In this context, as pointed out by Saltelli, issues such as relevance and transparency become critical [START_REF] Saltelli | Sensitivity analysis for importance assessment[END_REF]. Consequently, good modeling practice requires that the modeler provides an evaluation of the confidence in the model [START_REF] Crosetto | Sensitivity and uncertainty analysis in spatial modelling based on gis[END_REF]. This requires both uncertainty analysis in the model results and an evaluation of the part of output uncertainty due to each model input or a subgroup of the model inputs.

This work only focuses on the latter issue, which is the prerogative of sensitivity analysis (SA).

After quantifying uncertainties in a 0D arterial blood flow model of patient-specific upper body region including the circle of Willis in order to test its robustness and its relevance with respect to some choices and assumptions made in the model in previous works [START_REF] Lal | Data assimilation for identification of cardiovascular network characteristics[END_REF][START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF], this paper answers the question: which input variables or which input combinations are causing the largest effect on the model output? In other words, this work aims at addressing the following issues [START_REF] Pianosi | Sensitivity analysis of environmental models: A systematic review with practical workflow[END_REF]: which input parameters cause the largest variation in the output? Is there any input parameter whose variability has a negligible effect on the output? Are there interactions that amplify or dampen the variability induced by individual parameters? Specifically and in the words of Saltelli, this study shows how the uncertainty in the model output can be apportioned to different sources of uncertainty in the model inputs [START_REF] Saltelli | Sensitivity analysis for importance assessment[END_REF]. This type of approach, also known as importance measures or sensitivity indices, is known to be a variance-based method, which provides a factor-based decomposition of the output variance and implicitly assumes that this moment is sufficient to describe the output variability. Therefore, the main objective of our approach is to identify and rank in order of importance the important parameters in the blood flow model in order to establish a stopping criterion for the ensemble Kalman filter (EnKF) based algorithm used in previous works [START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF]. To be more precise, this work employs the Sobol' first order and total indices to analyze the link between the convergence of the 0D arterial blood flow model parameters estimated in previous works [START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF] using a data assimilation procedure and the importance of these parameters through SA with the following steps:

Proceed with the assimilation until the model with identified parameters fit the target observations, Identify the variables in the model which do not fully converge although target observations have been fitted in order to rank them by order of importance using the Sobol' sensitivity indices, Find a stopping criterion for the EnKF-based algorithm used in order to perform faster diagnosis in clinical cases.

Coupling EnKF and SA based on the Sobol' indices is not new. For instance, in Chen et al. [START_REF] Chen | Some sequential Monte Carlo techniques for Data Assimilation in a plant growth model[END_REF] SA based on the Sobol'indices is associated to EnKF for parameter screening purpose as part of a threestep prediction process. Also, nowadays it is widely recognized that the knowledge of sensitive input parameters can be beneficial for a range of purposes [START_REF] Pannell | Sensitivity analysis of normative economic models: theoretical framework and practical strategies[END_REF], such as better understanding of the relationships between input and output variables in a model, testing the robustness of a model results, development of robust and better model by reducing errors and uncertainties in the model, or by simplifying the model [START_REF] Hill | Practical use of computationally frugal model analysis methods[END_REF][START_REF] Lenhart | Comparison of two different approaches of sensitivity analysis[END_REF][START_REF] Bahremand | Distributed hydrological modeling and sensitivity analysis in torysa watershed, slovakia[END_REF][START_REF] Hill | Effective Groundwater Model Calibration: With Analysis of Data, Sensitivities, Predictions, and Uncertainty[END_REF]. Many publications on the subject, in particular the works of Saltelli et al. and references therein [START_REF] Saltelli | Sensitivity analysis for importance assessment[END_REF][START_REF] Saltelli | The role of sensitivity analysis in the corroboration of models and its links to model structural and parametric uncertainty (Guest editorial)[END_REF][START_REF] Saltelli | Sensitivity anaysis as an ingredient of modeling[END_REF][START_REF] Saltelli | On the relative importance of input factors in mathematical models[END_REF][START_REF] Saltelli | Global Sensitivity Analysis[END_REF], explain and illustrate these objectives.

The layout of the paper is as follows. The blood flow model under investigation, the biomedical problem to be addressed, and the detailed motivations for the current work are presented in Section 2. Section 3 is intended to introduce the sensitivity analysis method used in this work, namely sensitivity analysis based on the Sobol' indices. Sections 4 and 5 illustrate the usefulness and importance of the sensitivity indices in a patient-specific blood flow model, while Section 6 provides concluding remarks.

MATERIALS, PROBLEM DEFINITION, AND MOTIVATIONS

Patient-specific data and 0D blood flow model

The patient-specific data used in the current study comes mainly from medical imaging. These are essentially magnetic resonance angiography and magnetic resonance imaging (MRA&MRI) acquisitions, and a 3D time of flight magnetic resonance angiography (3D-TOF-MRA) images provided by the Department of Neuroradiology at the Centre Hospitalier Régional Universitaire de Montpellier (CHRU), Montpellier, France. Exploring these images with dedicated softwares, hemodynamic data, namely blood flow rates in ascending aorta (AA) and in the right and left internal carotid arteries (R-ICA & L-ICA), and morphological data including geometric measurements of some arteries, have been extracted. The missing geometric measurements were obtained from average data reported in the literature [START_REF] Alastruey | Modelling the circle of willis to assess the effects of anatomical variations and occlusions on cerebral flows[END_REF][START_REF] Reymond | Validation of a one-dimensional model of the systemic arterial tree[END_REF]. This also made it possible to adapt the arterial model proposed by Alastruey [START_REF] Alastruey | Modelling the circle of willis to assess the effects of anatomical variations and occlusions on cerebral flows[END_REF] to build the patient-specific arterial network as shown in Figure 1 (Figure 1A for the 1D model and Figure 1B for the equivalent 0D model). This patient-specific network consists of 33 arteries (each represented by one segment in Figure 1B) including the aorta, vertebral, carotid and brachial arteries together with the complete circle of Willis (coW) . For more details, the reader is referred to previous works [START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF].

Each segment in Figure 1B and in the 0D models in general is analogous to an electrical circuit [START_REF] Quarteroni | Coupling between lumped and distributed models for blood flow problems[END_REF][START_REF] Milišić | Analysis of lumped parameter models for blood flow simulations and their relation with 1d models[END_REF][START_REF] Pant | A methodological paradigm for patientspecific multi-scale cfd simulations: from clinical measurements to parameter estimates for individual analysis[END_REF][START_REF] Ursino | A model of cerebrovascular reactivity including the circle of willis and cortical anastomoses[END_REF], and is therefore represented as shown in Figure 2. For further details, please refer to previous studies [START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF].

Thus, applying Kirchhoff current and voltage laws to a single arterial segment assumed being filled with an incompressible Newtonian fluid, the following system of coupled ordinary differential equations for blood pressure and blood flow rates is obtained [START_REF] Quarteroni | Coupling between lumped and distributed models for blood flow problems[END_REF][START_REF] Milišić | Analysis of lumped parameter models for blood flow simulations and their relation with 1d models[END_REF][START_REF] Ursino | A model of cerebrovascular reactivity including the circle of willis and cortical anastomoses[END_REF] :

dP out dt = q in -q out C (Pressure equation) dq in dt = P in -P out -Rq in L (Flow equation), (1)
where for each arterial segment, P in , P out , q in , q out , R, L, and C are inlet blood pressure, outlet blood pressure, inlet flow rate, outlet flow rate, the viscous flow resistance, blood inertia, and the artery compliance, respectively. For arteries with a radius r < 0.2 cm, the inertial effect is ignored [START_REF] Olufsen | On deriving lumped models for blood flow and pressure in the systemic arteries[END_REF], and in this case, the flow equation is given by q in = (P in -P out )/R.

The electrical parameters (R, L, and C) of each arterial segment are related to the geometrical and physical parameters (r, h, l, and E) as follows [START_REF] Pant | A methodological paradigm for patientspecific multi-scale cfd simulations: from clinical measurements to parameter estimates for individual analysis[END_REF] : where r, h, l, µ, ρ, and E are the radius of the artery, the arterial wall thickness, the length of the arterial segment, the blood viscosity, the blood density, and the Young modulus, respectively.

At the bifurcations and at the junctions of the arterial trees shown in Figure 1, the enforcement of the mass and momentum conservation laws makes it possible to prescribe for diverging and merging flows the conditions at the nodes, namely that the blood pressure remains continuous and the blood flow rates must be conserved. For further details, the reader is referred to [START_REF] Sherwin | One-dimensional modelling of a vascular network in space-time variables[END_REF][START_REF] Gul | Mathematical modeling and sensitivity analysis of arterial anastomosis in the arm[END_REF][START_REF] Lal R | Data assimilation and uncertainty quantification in cardiovascular biomechanices[END_REF].

In order to take into account the effect of the downstream vasculature, the three-element Windkessel model (WK3-lumped parameter model) [START_REF] Reymond | Validation of a one-dimensional model of the systemic arterial tree[END_REF][START_REF] Boileau | A benchmark study of numerical schemes for one-dimensional arterial blood flow modelling[END_REF] is applied at the outlet of each terminal arterial segment.

As a reminder, the governing ordinary differential equation relating the instantaneous blood pressure and the blood flow rate for the WK3 model is written as follows:

dp(t) dt + p R D C = R P dq(t) dt + q(R P + R D ) R D C , (2)
where p, q, R P , R D , and C, are the instantaneous pressure at the inlet of the WK3 model, the instantaneous flow rate, the proximal resistance, the distal resistance, and the compliance of the vascular beds, respectively.

To solve the arterial network equations (the system of first-order ordinary differential equations consisting of equations ( 1) and ( 2)), the Fortran version of an implicit numerical integration solver DVODE [START_REF] Byrne | Stiff ode solvers: A review of current and coming attractions[END_REF][START_REF] Sp Corwin | Solving ODEs and DDEs with impulses[END_REF], available on http://www.radford.edu/ ˜thompson/vodef90web/, is used.

Problem definition

Given the physical model described above, this study focuses on the following optimization problem:

min x∈O ad J t, y(x,t), y obs (t) , (3)
x ∈ O ad ⊂ IR 21 , y(x,t) and y obs (t) ∈ IR 2 , where x and y obs (t) are independent variables and the time-dependent cost function J, of the form:

J t, y(x, t), y obs (t) = 1 2 y(x,t) -y obs (t) 2 .
Only the state variable, y(x,t), solution of a state equation F (y(x,t)) = 0 and the time variable, t, are related to the cost function J, the errors between the simulated and observed blow flow rates in the right internal carotid artery (R-ICA: compartment #21 in Figure 1A) and in the left internal carotid artery (L-ICA: compartment #23 in Figure 1A). The optimization parameter, x, which belongs to O ad , the optimization admissible domain [START_REF] Mohammadi | Applied Shape Optimization for Fluids (2nd Edition)[END_REF], is a vector of size 21 whose elements consist of 18 parameters (R Pi , R Di and C i where i = 9, 10, 19, 22, 29, 31 denotes the compartment numbers in Figure 1A) used to define the eleven Windkessel WK3 model boundary conditions in Figure 1 and

3 constants (k 1 , k 2 , k 3 )
defining the product of Young modulus, E, and thickness of arteries, h, by means of this empirical formula [START_REF] Olufsen | Structured tree outflow condition for blood flow in larger systemic arteries[END_REF]: Eh = r(k 1 e k2r + k 3 ). The observations variable, y obs (t), represents two blood flow rates time series obtained from MRA&MRI acquisitions at right and left internal carotid arteries as described Section 2.1. The state variable, y(x,t), is the estimation of the observation data, y obs (t), using the blood flow model presented above. The goal of this minimization problem is to find the right optimization parameter, x, which minimizes the error (the Euclidean distance) between the estimate of the observed data, y(x,t), and the observations themselves, y obs (t). More precisely, this consists in finding the right [START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF][START_REF] Saltelli | Global Sensitivity Analysis[END_REF][START_REF] Quarteroni | Coupling between lumped and distributed models for blood flow problems[END_REF][START_REF] Lal R | Data assimilation and uncertainty quantification in cardiovascular biomechanices[END_REF][START_REF] Byrne | Stiff ode solvers: A review of current and coming attractions[END_REF] that would allow to reproduce the flow rates observed in the patient-specific's right and left internal carotid arteries (ICAs).

set of 21 model parameters (k 1 , k 2 , k 3 , R Pi , R Di and C i where i = 9,
For this purpose, an ensemble Kalman filter based parameter estimation algorithm was developed and used. The latter was executed for 60 s (corresponding to 1,440 iterations, the maximum number of iterations arbitrary fixed) using an ensemble of size 30 and under assumptions. For more details either about these assumptions or about the parameter estimation algorithm and its execution, the reader is again referred to previous works [START_REF] Lal | Data assimilation for identification of cardiovascular network characteristics[END_REF][START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF].

Motivation: typical convergence using the EnKF-based parameter estimation algorithm

The data assimilation procedure based on the EnKF algorithm described and used in previous works [START_REF] Lal | Data assimilation for identification of cardiovascular network characteristics[END_REF][START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF][START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF] is an elegant way for parameter estimation. In these studies, the assimilation procedure was coupled either to a 1D arterial network or to the 0D compartment arterial network shown in Figure 1B with the aim of solving the optimization problem recalled above in (3). This made it possible to estimate hemodynamic parameters, which are difficult to identify noninvasively. For instance, in Lal et al. [START_REF] Lal | Data assimilation for identification of cardiovascular network characteristics[END_REF][START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF],

the Young modulus, E, and the optimization parameter, x, for the arterial network shown in Figure 1B were estimated as the solution of an inverse problem. This led, after assimilation, to patient-specific blood pressure assessment in cerebral arteries. How to explain, beyond the non-uniqueness of the optimization parameter detected in previous studies [START_REF] Lal | Data assimilation for identification of cardiovascular network characteristics[END_REF][START_REF] Lal | Non invasive blood flow features estimation in cerebral arteries from uncertain medical data[END_REF], the fact that the signals are retrieved whereas some variables have not fully converged? In other words, how to explain the poor identifiability of the optimization parameter?

Given also that iterating more does not change the signal overlap, when can the assimilation procedure be stopped to avoid unnecessary calculations? 0ne possible response (hypothesis) to the first issue raised above, namely the question of the poor identifiability of the optimization parameter -to be considered as an 'ansatz'-would be that the parameters not having fully converged, have little or no effect on the cost function, J, and therefore little or no effect in estimating the observed blood flow rate. One possible way to test this hypothesis is to resort to sensitivity analysis techniques in the hope of getting some clues about the influence of the estimated parameters on the observed blood flow rate estimate. The results obtained using this method could then be useful in addressing the second issue raised above. the variation in the output of a model can be attributed to variations of its input variables [START_REF] Pianosi | Sensitivity analysis of environmental models: A systematic review with practical workflow[END_REF]. The reader is referred to Saltelli for more details [START_REF] Saltelli | Global Sensitivity Analysis[END_REF][START_REF] Saltelli | Sensitivity Analysis[END_REF]. In general, SA is applied for multiple purposes, including support model calibration, verification and simplification, support the prioritization of efforts for uncertainty reduction, or help with model-based decision-making [START_REF] Pianosi | Sensitivity analysis of environmental models: A systematic review with practical workflow[END_REF][START_REF] Norton | An introduction to sensitivity assessment of simulation models[END_REF][START_REF] Pianosi | A matlab toolbox for global sensitivity analysis[END_REF][START_REF] Song | Global sensitivity analysis in hydrological modeling: Review of concepts, methods, theoretical framework, and applications[END_REF]. Such purposes are generally implemented as four main objectives of SA [START_REF] Saltelli | Global Sensitivity Analysis[END_REF][START_REF] Iooss | Introduction to sensitivity analysis. Handbook of uncertainty quantification[END_REF]: the Factors Fixing (FF) setting or screening [START_REF] Pianosi | Sensitivity analysis of environmental models: A systematic review with practical workflow[END_REF][START_REF] Sarrazin | Global sensitivity analysis of environmental models[END_REF], the Factors Prioritization (FP) setting or ranking [START_REF] Van Werkhoven | Characterization of watershed model behavior across a hydroclimatic gradient[END_REF][START_REF] Sin | Global sensitivity analysis in wastewater treatment plant model applications: Prioritizing sources of uncertainty[END_REF][START_REF] Hartmann | Testing the realism of model structures to identify karst system processes using water quality and quantity signatures[END_REF], the Variance Cutting (VC) setting [START_REF] Saltelli | On the relative importance of input factors in mathematical models[END_REF][START_REF] Iooss | Introduction to sensitivity analysis. Handbook of uncertainty quantification[END_REF], and the Factors Mapping (FM) [START_REF] Iooss | Introduction to sensitivity analysis. Handbook of uncertainty quantification[END_REF][START_REF] Spear | Eutrophication in peel inletii. identification of critical uncertainties via generalized sensitivity analysis[END_REF][START_REF] Singh | A vulnerability driven approach to identify adverse climate and land use change combinations for critical hydrologic indicator thresholds: Application to a watershed in pennsylvania, usa[END_REF]. These SA objectives and methods in turn can be gathered into two main categories [START_REF] Iooss | Introduction to sensitivity analysis. Handbook of uncertainty quantification[END_REF]: Local Sensitivity Analysis (LSA), which is interested in calculating some kind of derivative or partial derivatives (gradient) of the model at a specific point of the input variable space [START_REF] Iooss | Introduction to sensitivity analysis. Handbook of uncertainty quantification[END_REF][START_REF] Turányi | Sensitivity analysis of complex kinetic systems. tools and applications[END_REF][START_REF] Turányi | Sensitivity analysis of complex kinetic systems. tools and applications[END_REF][START_REF] Hill | Evaluating Estimated Parameter Values and Parameter Uncertainty[END_REF][START_REF] Ljung | System identification[END_REF][START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF][START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF], and Global Sensitivity Analysis (GSA), which unlike LSA rather attempts to explore the entire space of the input factors. [START_REF] Saltelli | Sensitivity anaysis as an ingredient of modeling[END_REF][START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF][START_REF] Zhou | Global sensitivity analysis[END_REF]. The Sobol' indices that will be used in this work belong to this category.

The Sobol' sensitivity indices

The literature on sensitivity measures is abundant [START_REF] Saltelli | Sensitivity analysis for importance assessment[END_REF][START_REF] Saltelli | Sensitivity Analysis[END_REF][START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF][START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF][START_REF] Archer | Sensitivity measures,anova-like techniques and the use of bootstrap[END_REF][START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their monte carlo estimates[END_REF][START_REF] Sobol | Derivative based global sensitivity measures and their link with global sensitivity indices[END_REF][START_REF] Saltelli | Sensitivity analysis practices: Strategies for modelbased inference[END_REF][START_REF] Melis | Bayesian sensitivity analysis of a 1d vascular model with gaussian process emulators[END_REF]. According to Sudret, among this vast literature on sensitivity measures, the Sobol' sensitivity indices have received much attention since they provide accurate information for most models [START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF]. Indeed, first introduced by Sobol in 1990 [START_REF] Sobol | On sensitivity estimation for nonlinear mathematical models[END_REF] and based on his earlier work on the Fourier Haar series [START_REF] Sobol | Multidimensional quadrature formulas and haar functions. Izdat[END_REF], the indices were developed for the purpose of estimating not only the functional structure of Y (the model output) but also the sensitivity of Y with respect to any model input variable or any subgroup of model input variables. Other techniques existing in the literature, such as active subspaces methods (see Constantine et al. [START_REF] Constantine | Active subspace methods in theory and practice: applications to kriging surfaces[END_REF], Palar et al. [START_REF] Palar | Exploiting active subspaces in global optimization: How complex is your problem[END_REF], and references therein for more details), share the same goal as the Sobol' indices based approach and could also have been used. However in this Section, only the essential features of the Sobol' sensitivity indices method as presented by Sobol [START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF] and Archer et al. [START_REF] Archer | Sensitivity measures,anova-like techniques and the use of bootstrap[END_REF] will be followed.

Expansion into summands of increasing dimensions

Let I the unit interval [0, 1], I d = {X : 0 ≤ X i ≤ 1, i = 1, . . . ,
f : I d -→ IR p X -→ Y = f (X) = f (X 1 , . . . , X d ).
According to Sobol, the function f can be uniquely expanded in a series of 2 d summands of increasing dimensions [START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF][START_REF] Archer | Sensitivity measures,anova-like techniques and the use of bootstrap[END_REF][START_REF] Sobol | On sensitivity estimation for nonlinear mathematical models[END_REF][START_REF] Hoeffding | A class of statistics with asymptotically normal distribution[END_REF]:

f (X) = f 0 + d i=1 f i (X i ) + 1≤i<j≤d f ij (X i , X j ) + • • • + f 12...d (X 1 , . . . , X d ). (4)
provided that f 0 is constant and

I f i1•••is (X i1 , . . . , X is ) dX i k = 0, 1 ≤ k ≤ s.

Functional decomposition of variance (ANOVA-representation) and Sobol' indices estimation

The total variance D of f (X) can be obtained and decomposed in the same way as Eq.( 4) [START_REF] Archer | Sensitivity measures,anova-like techniques and the use of bootstrap[END_REF][START_REF] Sobol | On sensitivity estimation for nonlinear mathematical models[END_REF][START_REF] Efron | The jackknife estimate of variance[END_REF]:

D = d i=1 D i + 1≤i<j≤d D ij + • • • + D 12...d ,
where D i1...is =

I d f 2 i1...is (X i1 , • • • , X is ) dX i1 .
. . dX is is the contribution to the total variance from term f i1...is in Eq.( 4). The so-called Sobol' indices denoted by S i1...is are the fractions of the total variance of f (X) due to any individual summand in Eq (4). They are obtained as follows [START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF]:

S i1...is = Di 1 ...is D
and their sum over all possible combinations of indices is 1. The number of sensitivity indices thus calculated, from first order to order d, is 2 d -1. Therefore, when the number of input variables d becomes too large, the number of sensitivity indices become exceedingly large so that estimating and interpreting all these indices quickly becomes impossible. That is why the practitioner does not usually estimate indices of order higher than two. However, another alternative was introduced by Homma and Saltelli [START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF] with the total sensitivity indices.

The total sensitivity indices

Also called "total indices" or "total effects", these indices express the total sensitivity of the variance of Y (the model output) with respect to any given input variable, i.e. the sensitivity with respect to this variable in all its forms (sensitivity with respect to the variable alone and sensitivity with respect to the interactions of this variable with other variables). The total sensitivity index, S Ti , with respect to the variable X i is therefore defined as the sum of all sensitivity indices in regard to the variable X i and can be written as follows [START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF][START_REF] Iooss | A Review on Global Sensitivity Analysis Methods[END_REF]:

S Ti = l∈@i S l
where @i represent all the admissible subsets of {1, • • • , d} including i. For instance, for a model with three input variables, the total index with respect to the variable X 1 is S T1 = S 1 + S 12 + S 13 + S 123 .

Numerical estimation of the Sobol'sensitivity indices

For numerical computation of the Sobol' sensitivity indices, several techniques and algorithms have been developed and can be found in the literature. For instance, the reader is referred to Sobol [START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF][START_REF] Sobol | On sensitivity estimation for nonlinear mathematical models[END_REF] and Saltelli [START_REF] Saltelli | Making best use of model evaluations to compute sensitivity indices[END_REF] where Monte Carlo sampling based methods are proposed either for first order and interaction indices computation [START_REF] Sobol | Sensitivity analysis for nonlinear mathematical models[END_REF] or for first order and total indices computation [START_REF] Saltelli | Making best use of model evaluations to compute sensitivity indices[END_REF]. Given the cost of such methods in terms of number of model calls (the rate of convergence is in O(N -1/2 ) where N represents the sample size) [START_REF] Iooss | A Review on Global Sensitivity Analysis Methods[END_REF], many extensions of these methods and other alternative methods for reducing the complexity of indices computation have been proposed by several researchers [START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF][START_REF] Mckay | Comparison of three methods for selecting values of input variables in the analysis of output from a computer code[END_REF][START_REF] Sobo | Quasi-monte carlo methods[END_REF][START_REF] Niederreiter | Random Number Generation and Quasi-Monte Carlo Methods[END_REF][START_REF] Owen | Randomly permuted (t,m,s)-nets and (t, s)-sequences. Monte Carlo and Quasi-Monte Carlo Methods in Scientific Computing[END_REF][START_REF] Owen | Monte carlo extension of quasi-monte carlo[END_REF][START_REF] Owen | Monte carlo, quasi-monte carlo, and randomized quasi-monte carlo[END_REF][START_REF] Saltelli | A quantitative model-independent method for global sensitivity analysis of model output[END_REF][START_REF] Saltelli | How to avoid a perfunctory sensitivity analysis[END_REF][START_REF] Cukier | Nonlinear sensitivity analysis of multiparameter model systems[END_REF][START_REF] Tissot | Bias correction for the estimation of sensitivity indices based on random balance designs[END_REF][START_REF] Tarantola | Estimating the approximation error when fixing unessential factors in global sensitivity analysis[END_REF][START_REF] Janon | Asymptotic normality and efficiency of two sobol index estimators[END_REF][START_REF] Veiga | Local polynomial estimation for sensitivity analysis on models with correlated inputs[END_REF].

These proposed extensions include pseudo-probabilistic simulation methods such as Latin Hypercube Sampling (LHS), Fourier Amplitude Sensitivity Test (FAST) based methods, Local polynomial and polynomial chaos expansions based methods, Quasi-Monte Carlo (QMC) based methods, etc., and their variants. The details of these different methods are beyond the scope of this paper. Indeed, in the present work, the package sensitivity of R environment available on the CRAN website (https://cran. r-project.org/web/packages/sensitivity/sensitivity.pdf) [START_REF] Iooss | Package sensitivity, version 1.17.1. The Comprehensive R Archive Network[END_REF] will be used to compute the Sobol' sensitivity indices. The objective will be to estimate both first order and total indices in order to be able to take into account the interactions of each parameter with all the others. Thus, the function sobolmartinez of the package sensitivity, which implements the Monte Carlo estimation of the Sobol' indices for both first-order and total indices, will be used.

APPLICATION TO THE PATIENT-SPECIFIC 0D BLOOD FLOW MODEL

The performance of the Sobol' indices is investigated using the blood flow model described in Section 2 with the aim of confirming or refuting our initial hypothesis formulated in Section 2.3 about the poor identifiability of the optimization parameter. But before that, let us recall the procedure for generating the patient-specific database.

Database generation

During the EnKF-based algorithm execution, as early presented in Section 2.3, one noticed that from 8.35 s, there was a very good agreement between the observed and predicted blood flow waveforms and also that iterating more with the assimilation algorithm did not modify the recovery of the observed blood flow rates (see Figure 4). This means that after 8.35 s (about 10 cardiac cycles, or about 200 EnKF iterations), the right set of model parameters allowing to reproduce the blood flow rates observed in the ICAs is found and therefore the execution of the algorithm can be stopped. But, to ensure not to stop the algorithm too early, we let it run until the end (i.e. up to the 1440th iteration) and only saved the results for the last 20 EnKF iterations, namely the optimization parameter, x (the 21 estimated variables with their associated uncertainties), the cost function, J, and also the variability of the cost function with respect to each estimated variable at each EnKF-based algorithm iteration. These are used to make up the patient-specific database and to build a simplified model of the patient-specific blood flow model in order to be able to easily and properly perform SA.

Sensitivity analysis and hypothesis testing using the Sobol' indices

To identify both the most influential, the weakly influential, and the non-influential input parameters in the patient-specific blood flow model in order to address the question of the reliability of the original assumption stated above in Section 2.3, SA using the Sobol' indices is performed on a simplified model of the patient-specific blood flow model. This is to simplify the Sobol' indices calculations using simple information from final iterations of the EnKF procedure. Our sensitivity analysis is therefore under the hypothesis of local linearity as this analysis is performed after convergence and needs to be valid around the final solution, and not over the whole variable space. It is classical to locally substitute nonlinear models with their linearized versions for sensitivity and/or stability analysis, through spectral analysis for instance. Here, the substitution permits to keep the Sobol' indices calculations simple and possible a posteriori. Thus, the 21 variables making up the optimization parameter are assumed to be independent and their relationship with the cost function, J, is assumed to be governed by the this simple linear model based on the saved data in the last 20 EnKF iterations and defined for each of these iterations as follows:

L : IR 21 -→ IR 2 X -→ L(X) = BX = J = (J RICA , J LICA ) ,
where X, L, and J stand for the optimization parameter (also referred as x throughout the paper), the linear mapping associated to the matrix B, and the errors between the simulated and observed blow flow rates in R-ICA and L-ICA, respectively, in each of the last 20 EnKF iterations. The components of J, J RICA and J LICA , represent the cost function obtained when using only the blood flow rate in R-ICA as observed data and the cost function obtained when using only the blood flow rate in L-ICA as observed data, respectively, in each of the last 20 EnKF iterations. For each of these iterations, B is a 2-by-21 matrix obtained by assembling the local variability values of the components of the cost function, J, with respect to each component of the optimization parameter. This means that for each iteration considered,

B = ∂J RICA ∂X i , ∂J LICA ∂X i , i = 1, • • • , 21,
where ∂J RICA

RESULTS AND STOPPING CRITERION SETTING

Results and discussion

Figure 5 shows typical results of the first order and total indices over some of the last 20 EnKF iterations when using only either J RICA or J LICA as the simplified model output. It is important to note that similar results are obtained when considering the last 40, 60, 80, 100, 120 iterations and even beyond (the only constraint is not to go as far as considering the first 200 iterations, which provide rather disparate results).

The interpretation of these results is not entirely straightforward, obvious. ) with R P31 in the lead. This means that R P31 seems to be the most influential parameter, R P10 , R P19 , R P22 , and R P29 , the slightly influential ones and the rest seems to be non-influential. In order to review this general trend and to provide for each estimated parameter single representative main and total indices, which are much more easy to interpret, the means and standard deviations of the first order and total indices over all the last 20 EnKF iterations for each of the 21 estimated parameter have been calculated. The results are shown in Figure 6 when using either only J LICA (Figure 6A) or only J RICA (Figure 6B) as the simplified model output. The results from these Figures are perfectly consistent with the general trend identified above and therefore confirm it. Indeed, from the main and total mean indices reported in Figure 6, five indices measuring the importance of the five proximal resistances cited above as influential parameters stand out from all others with the same order of importance. Comparing these results with those illustrated in Figure 3, it can be seen that the parameters found not to be fully convergent in this Figure have been identified here as non-influent parameters. Furthermore, it is also noted that all the parameters identified here as influential had fully converged, as shown in Figure 3. This is consistent with the hypothesis formulated from the outset and therefore, confirms and validates it at the same time.

However, it is still legitimate to ask the question of the reliability, the relevance or the credibility of these results. To ensure that these results about the importance of estimated parameters make sense, the factors fixing setting and the factors prioritization setting are performed. In concrete terms, the EnKFbased parameter estimation algorithm is executed using the following two scenarii with regard to the set of model input parameters:

Scenario A: The algorithm is executed with a set of input parameters created by keeping unchanged the default initialization values of the five important parameters detected above, and fixing during all EnKF iterations all other parameters to their nominal values taken in previous works [START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF].

Scenario B: The algorithm is now executed using a set of input parameters created by fixing during all EnKF iterations the five important parameters detected above to their nominal values (estimated values) taken in previous works [START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF] and keeping unchanged the default initialization values of all other parameters.

At the end, for each of both scenarii considered, the normalized cost function over time, J norm , which is obtained by dividing the cost function at each EnKF iteration, J by iter., by the cost function at the first EnKF iteration, J of 1 st iter., is compared to the reference normalized cost function i.e., that previously obtained using the set of the default input parameters [START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF]. The results of this comparison are shown in This corroborates both the fact that the input parameters fixed in scenario A are really unessential in the blood flow rates estimation in the patient-specific ICAs and the fact that only the input parameters fixed in scenario B are important in the blood flow rates estimation in the patient-specific ICAs. This means that for the purpose of the blood flow rate estimation in this patient-specific ICAs, the blood flow model can be reduced to the five input parameters detected as influential instead of twenty-one and still well performs.

Choice of a stopping criterion in the assimilation procedure

With the dust settled over the parameter importance issue, the essential question remains as to how a stopping criterion can be fixed for the EnKF-based parameter estimation algorithm, which under the current situation needs 7 hours on a 3GHz CPU for one inference, in order to prevent from unnecessary over-solving. Many possibilities can be considered. However, in view of the sensitivity and convergence analysis performed above, it seems judicious to propose, for the blood flow model examined here, a stopping criterion based only on the convergence of the most important parameters. Thus, the following stopping criterion is proposed: the EnKF-based algorithm execution is automatically interrupted as soon as, at a given iteration, all the absolute values of the differences between successive values of the five important parameters each reach the threshold value of 10 -4 . The choice of this threshold value is based on the currently available values of the difference between consecutive values of the estimated parameters after 8.35 s (about 200 EnKF iterations). This typically leads to a stopping time t = 10 s, reducing by a factor of six the time to solution. This means that in most situations the calculation (one inference on a 3GHz CPU) is now performed in about an hour. Therefore, the main advantage of establishing this stopping criterion is that it can be used to achieve faster diagnosis in clinical situations. However, this choice of the stopping criterion is too much adapted to the particular situation of interest here. Indeed, it shows how the current sensitivity analysis permits to decide when to stop based on how the most important variables behave: here the resistances but following the outcome of our sensitivity analysis, for another patient, it could be other variables. Nevertheless, several tests not illustrated here for the sake of brevity show that the results obtained are quite robust with respect to the choice of the last 20 iterations for sensitivity analysis. Indeed, these tests show that the results of the Sobol' indices computation were invariant (i.e., similar to the case where only the last 20 EnKF iterations have been considered) when considering the last 40, 60, 60, 100, 120 EnKF iterations and even much more, whereas they are very different and disparate before the first 200 iterations. This robustness actually somehow shows that the proposed stopping criterion is not aberrant. In any case, the integration of this stopping criterion to the current algorithm allow to have a parameter estimation procedure much more efficient for future tests.

R = 8µl πr 4 (

 4 Hagen-Poiseuille equation), L = ρl πr 2 , and C = 3πr 3 l 2Eh ,

4 Figure 1 .

 41 Figure 1. The network of a 1D blood flow model of the upper body arteries and of the CoW (A), and its equivalent compartment model (B), adapted from [20]. The lines size and the arrows indicate the relative size of the arteries and the direction of flow, respectively. The eleven red circles in each figure represent Windkessel WK3 model boundary conditions.

Figure 2 .

 2 Figure 2. A single compartment circuit illustration.

  Figures 3E and 3F, respectively), from 8.35 s (i.e., after 10 cardiac cycles), there is a very good agreement between the target and predicted flow rate waveforms meaning that the observations are well recovered (Figure3A,B). It is also noted that iterating more with the assimilation algorithm does not modify the recovery of the signals (Figure4C,D), so we would like to stop the assimilation procedure earlier to avoid unnecessary calculations. These two findings raise the following interesting issues, which are worth exploring, examining, more closely:

Figure 3 .Figure 4 .

 34 Figure 3. Typical convergence of some parameters previously estimated in [10] using the ensemble Kalman filter (EnKF) based algorithm. A,B: Comparison of the model simulated blood flow rate waveform in R-ICA (A) and L-ICA (B) to the target signals. C: Mean and standard deviation of estimated Young moduli for three of the 33 segments of the arterial networks in Figure 1. The evolution is over 8.35 s corresponding to about 10 cardiac cycles. D,E,F : Histories of the two first moments of some estimated parameters (Windkessel WK3 model parameters in Figure 1) during the EnKF-based algorithm iterations.

  ∂Xi and ∂J LICA ∂Xi refer to the local variability values of the cost function, J RICA , and to the local variability values of the cost function, J LICA , respectively, with respect to each component of the optimization parameter. The Sobol' first order indices and total indices for each of the 21 variables making up the optimization parameter in each of the last 20 EnKF iterations are estimated using the function sobolmartinez from the R package sensitivity with the following arguments: the simplified model defined above and two 10,000-random samples generated following a multivariate normal distribution with means equal to mean values of the 21 estimated parameters at each EnKF iteration considered and covariance matrix whose the elements consist of the standard deviations associated to the 21 estimated parameters in each EnKF iteration. This assumption of normal distribution for the input parameters is intended to produce samples compatible with the fundamental Gaussian hypothesis used in the EnKF procedure [8-10].

Figure 7 .

 7 Figure 7. One can notice that there is virtually no significant change on the reference normalized cost function J norm in scenario A. However, a reduction of the model output uncertainty can be observed in scenario B. Indeed, the values of the normalized cost function over time in scenario A and the reference ones are almost identical while those obtained in scenario B are relatively lower than the reference ones.

Figure 5 .

 5 Figure 5. First order and total indices for the 21 estimated parameters over some of the last 20 EnKF iterations saved and using either only the cost function, J LICA (left), or only the cost function, J RICA (right), as the simplified model output. This is a representative sample of typical Sobol' indices estimation results for the 21 parameters obtained over the last 20 iterations. The indices obtained over other iterations, which are not presented here, show similar patterns to those illustrated in this figure.

Figure 6 .

 6 Figure 6. Means and Standard deviations over the last 20 EnKF iterations of the first order and total indices for the 21 estimated parameters using at each EnKF iteration either only the cost function, J LICA (A), or only the cost function, J RICA (B), as the simplified model output.

Figure 7 .

 7 Figure 7. Comparison between the normalized cost function, Jnorm, obtained in[START_REF] Rapadamnaba | Backward sensitivity analysis and reducedorder covariance estimation in noninvasive parameter identification for cerebral arteries[END_REF] (referred to as Previous work) with those obtained in both scenarii considered: (A) for scenario A and (B) for scenario B. The solid coloured line shows the normalized cost function value over time (i.e., the cost function, J, over all ensemble Kalman filter (EnKF) algorithm iterations divided by the cost function, J, at the first EnKF algorithm iteration).

  It is even very delicate and should therefore be done with care. Indeed, from one Figure to another, diverse and various patterns are visible. However, looking very carefully and closely at these Figures, a general trend emerges. The blood flow model behaves as if it depends upon only a few parameters, with the rest having comparatively little influence. More precisely, only five of the 21 estimated parameters appear to be important. These are all proximal resistance parameters (R P10 , R P19 , R P22 , R P29 , R P31

R. RAPADAMNABA ET AL.

ACKNOWLEDGMENT

The authors gratefully acknowledge Prof.-Dr. V. Costalat, Dr. E. Le Bars and Dr. J. Deverdun from the Department of Neuroradiology of the Centre Hospitalier Régional Universitaire de Montpellier, Gui de Chauliac, Montpellier, France and Dr. J. Siguenza for providing the medical data.