Sébastien Roux
email: sebastien.roux@inrae.fr

Patrice Loisel

Samuel Buis

Inverting Regional Sensitivity Analysis to reveal sensitive model behaviors

Keywords: Multivariate sensitivity analysis, Target sensitivity analysis, Sobol' indices, Cluster-based GSA, Graphical sensitivity analysis

We address the question of sensitivity analysis for model outputs of any dimension using Regional Sensitivity Analysis (RSA). Classical RSA computes sensitivity indices related to the impact of model inputs variations on the occurrence of a target region of the model output space. In this work, we invert this perspective by proposing to find, for a given target model input, the region whose occurrence is best explained by the variations of this input. When it exists, this region can be seen as a model behavior which is particularly sensitive to the variations of the model input under study. We name this method iRSA (for inverse RSA). iRSA is formalized as an optimization problem using region-based sensitivity indices and solved using dedicated numerical algorithms. Using analytical and numerical examples, including an environmental model producing time series, we show that iRSA can provide a new graphical and interpretable characterization of sensitivity for model outputs of various dimensions.

Introduction

The analysis of models with complex outputs (temporal, spatial, heterogeneous) is one of the current challenges of sensitivity analysis [START_REF] Marrel | Sensitivity analysis of spatial and/or temporal phenomena[END_REF][START_REF] Saltelli | Sensitivity analysis: A discipline coming of age[END_REF]. There are two main approaches to handle multidimensional outputs. The first one aims at computing a set of scalar indices associated to coordinates along a basis of functions on which the model outputs are projected ( [START_REF] Campbell | Sensitivity analysis when model outputs are functions[END_REF][START_REF] Lamboni | Multivariate global sensitivity analysis for dynamic crop models[END_REF][START_REF] Xiao | Multivariate global sensitivity analysis for dynamic models based on wavelet analysis[END_REF]). Its flexibility and interpretability are limited by the choice of the basis. The second approach summarizes the impact of the model inputs on the variability of all the model outputs using a single index [START_REF] Lamboni | Multivariate sensitivity analysis to measure global contribution of input factors in dynamic models[END_REF][START_REF] Gamboa | Sensitivity analysis for multidimensional and functional outputs[END_REF][START_REF] Lamboni | Multivariate sensitivity analysis: Minimum variance unbiased estimators of the firstorder and total-effect covariance matrices[END_REF][START_REF] Xu | Generalized sensitivity indices based on vector projection for multivariate output[END_REF]. It provides useful aggregated indices but does not allow a fine understanding of model outputs sensitivity.

Spear and Hornberger introduced in [START_REF] Spear | Eutrophication in peel inlet-ii. identification of critical uncertainties via generalized sensitivity analysis[END_REF] the concept of model behavior in sensitivity analysis through the Regional Sensitivity Analysis (RSA) approach. The principle is to start from the definition of a target region of the output space (denoted as "behavioral") and to analyze the impact of the variations of model inputs on its occurrence. Using model behaviors expressed as regions of the output space appears to be an efficient method to get interpretable characterizations of model properties. It has also the property to scale to any dimension of the output space [START_REF] Roux | Cluster-based gsa: Global sensitivity analysis of models with temporal or spatial outputs using clustering[END_REF]. Among the last developments on RSA, two are of particular interest in the present study: i) the application of RSA in the context of reliability engineering to characterize parameter sensitivity in relation to a critical domain of the output space (e.g. the failure domain of a system) using sensitivity measures compatible with rare events and taking into account interactions (Target SA, [START_REF] Marrel | Statistical developments for target and conditional sensitivity analysis: application on safety studies for nuclear reactor[END_REF][START_REF] Idrissi | Developments and applications of shapley effects to reliabilityoriented sensitivity analysis with correlated inputs[END_REF]), ii) its application in combination with a clustering procedure in order to characterize parameter sensitivity with respect to the dominant model behaviors detected in the model output space (Distance-Based Generalized SA [START_REF] Fenwick | Quantifying asymmetric parameter interactions in sensitivity analysis: application to reservoir modeling[END_REF], Cluster-based GSA, [START_REF] Roux | Cluster-based gsa: Global sensitivity analysis of models with temporal or spatial outputs using clustering[END_REF]).

These approaches rely on an a priori characterization of the behaviors (regions of the output space) to be analyzed. They are identified either by experts or automatic clustering of the simulations. In this work, we propose a new perspective on the link between behaviors and sensitivity analysis based on an extension of the above mentioned developments on RSA. Instead of trying to a priori characterize target regions of the output space, we propose to use an optimization procedure in order to reveal the regions of the output space the most sensitive to the variations of a given input, i.e. whose occurrences are the best explained by the variations of this input. We name this approach iRSA (for inverse Regional Sensitivity Analysis). The formalization of iRSA in terms of principle and numerical algorithms will be presented in Section 2. The results of the method application are presented in Section 3 on three examples: a model with scalar outputs which is analytically solved and two models on which numerical algorithms are tested, a model with 2D outputs and an environmental model producing time series.

iRSA methods

In this section, we formalize iRSA in a general context, with a focus on the sensitivitybased optimization criteria that can be used and on their minimization in the context of multivariate outputs.

Notations

We consider a model f whose inputs are noted X i , with i ∈ [1, n] and have independent distributions. The model output is noted Y = f (X 1 , .., X n ) and can be multivariate.

Model behaviors in the context of RSA are defined as regions of the model output space. Depending on the formalization context, these regions are continuous sets (as we will see in Section 3.1) or discrete set (in Sections 3.2 and 3.3). In the latter case, they are ensembles of simulated points and we will denote them as clusters of model outputs.

In the following, we generically denote by C the model output space.

iRSA principle

The principle of iRSA is to find for each input factor X i a partition of the model output space C that maximizes a cost function related to parameter sensitivity. We will be interested in either finding the region C whose occurrence is the most influenced by a given input, in which case iRSA will provide the partition (C, C) (where C is the complement of C), or the two non-overlapping regions C and C whose transition from the one to the other is the most influenced by a given input, in which case iRSA will provide the partition made of three regions (C, C , C ∪ C ).

Main criteria and associated interpretations

We consider the region-based sensitivity indices introduced in [START_REF] Roux | Cluster-based gsa: Global sensitivity analysis of models with temporal or spatial outputs using clustering[END_REF]. These indices use membership functions which characterize the level of membership of any point of the output space to a given region C. In this work, we consider binary membership functions that can be expressed using the indicator function 1 C (.).

The idea behind these indices is simply to build sensitivity indices from the standard Sobol' indices [START_REF] Sobol | Sensitivity estimates for nonlinear mathematical models[END_REF][START_REF] Saltelli | Sensitivity analysis[END_REF] obtained when transforming the multivariate output into a scalar one using the membership functions.

Indices based on single region membership

The first and simplest region-based indices that can be used for an input X i are the Sobol' indices (first or total: SI C i and TSI C i ) associated to the membership function of a region C defined over the model output space C.

SI C i = V [E [1 C (Y )|X i ]] V [1 C (Y )] TSI C i = 1 - V [E [1 C (Y )|X ∼i ]] V [1 C (Y )]
In that case, iRSA aims at finding for each input X i a partition (C * i , C * i ) so that C * i maximizes either the first or the total region-based Sobol' index of input X i . By definition, the partition (C * i , C * i ) defines the two behaviors whose transition from the one to the other is the most influenced by the model input X i . C * i (or equivalently C * i ) can also be seen as the region whose occurrence is best explained by the variations of this input.

The optimization problems can be written :

C * i = arg max C∈C SI C i (SI C i based iRSA) C * i = arg max C∈C TSI C i (TSI C i based iRSA) 3
2.3.2. Indices based on cluster membership differences A extension of these indices was also proposed in [START_REF] Roux | Cluster-based gsa: Global sensitivity analysis of models with temporal or spatial outputs using clustering[END_REF] by considering Sobol' indices on differences of membership functions.

SI CC i = V [E [1 C (Y ) -1 C (Y )|X j ]] V [1 C (Y ) -1 C (Y )] TSI CC i = 1 - V [E [1 C (Y ) -1 C (Y )|X ∼i ]] V [1 C (Y ) -1 C (Y )]
In that case, iRSA aims at finding for each input The optimization problems can be written :

X i a partition (C * i , C * i , C * i ∪ C * i ) so that that (C * i , C * i )
(C 1 * i , C 2 * i ) = arg max C,C ∈C,C∩C =∅ SI CC i (SI CC i based iRSA) (C 1 * i , C 2 * i ) = arg max C,C ∈C,C∩C =∅ TSI CC i (TSI CC i based iRSA)

Generic optimization algorithms

The previous criteria can be optimized analytically only for simple models (an example is presented in Section 3.1). In a general setting, numerical approaches must be considered. For such numerical algorithms, the output space C is the set of all simulated points. Regions of C can thus be seen as clusters of simulated points and iRSA as a clustering algorithm that uses sensitivity-based criteria.

All proposed algorithms require that a classical sensitivity workflow compatible with scalar outputs can be implemented on the model of interest and, more precisely, that a design of experiment X is available for the sensitivity analysis. We note Y the matrix of multivariate outputs obtained when applying the model on the design matrix X. We suppose that applying the model on X allows to compute Sobol' type sensitivity indices on scalar outputs (for example, X is defined using a Jansen pick and freeze method that allows the computation of Sobol' indices using the Jansen estimation formula [START_REF] Jansen | Analysis of variance designs for model output[END_REF]).

In the following, we propose a set of algorithms implementing iRSA. Their shared principle is to perform a pre-clustering of C into K Y clusters noted (C 1 , ..., C K Y ) in order to reduce the computational cost. This pre-clustering is performed using a classical distance-based approach (denoted as ClustF un Y in the different algorithms), typically a K-Means algorithm or a hierarchical clustering method, but other methods can be considered to better suit a given application.

Algorithm for indices based on single cluster membership

We present a generic algorithm working for any region-based sensitivity criterion Crit provided that it is defined using single membership functions. The most simple criteria are SI C i or TSI C i as presented in the previous sections. The algorithm has several parameters: the index i of the input under study, the design matrix X and the associated output matrix Y , the criterion on single membership function Crit, the pre-clustering method ClustF un Y , the number K Y of clusters for the pre-clustering step and a size parameter γ that can be used to avoid partition made of too small regions. The optimization of Crit is performed by enumerating all partitions of (C 1 , .., C K Y ) into two non-empty sets. The total number of such 2-partitions is the Stirling number S K Y 2 = 2 K Y -1 -1 which drives the complexity of the algorithm. The algorithm returns the best partition under the size constraint as well as the associated sensitivity score. We denote as iRSA SM (SM for for Single Membership) this algorithm:

Algorithm iRSA SM(i, X, Y , Crit, ClustF un Y , K Y , γ) Apply ClustF un Y (K Y ) on Y to get K Y clusters (C 1 , .., C K Y ) for all 2-partitions (P k , P k ) k∈[1:S K Y 2 ] of (C 1 , .., , C K Y ) do Compute γ k = min Card(P k ), Card(P k ) Compute Crit(i, X, Y , P k ) end for Get k * = arg max k∈[1:S K Y 2 ],γ k ≥γ Crit(i, X, Y , P k ) return (P k * , P k * ) and Crit(i, X, Y , P k * )
Note that by storing all sensitivity indices for each partition in the loop, it is possible to process all inputs at once during the post-processing step, thus avoiding to run the loop again for different inputs.

Algorithm for indices based on cluster membership differences

In the same spirit as the previous algorithm, the optimization of a sensitivity criterion based on membership function differences such as SI CC i or TSI CC i (and denoted generically as Crit in the detail of the algorithm), is performed by enumerating all partitions. Due to the structure of these indices, it is necessary to enumerate all partitions of (C 1 , .., C K Y ) into three non-empty sets. The total number of such 3-partitions is the Stirling number

S K Y 3 = 1 6 3 K Y -3 • 2 K Y + 3 .
We denote as iRSA DM (DM for Difference of Membership) this algorithm.

Algorithm iRSA DM(i, X, Y , Crit, ClustF un Y , K Y , γ) Apply ClustF un Y (K Y ) on Y to get K Y clusters (C 1 , .., C K Y ) for all 3-partitions Q k 1 , Q k 2 , Q k 3 k∈[1:S K Y 3 ] of (C 1 , .., , C K Y ) do Compute γ k = min Card(Q k 1 ), Card(Q k 2 ), Card(Q k 3 ) Compute for (p, q) ∈ I = {(1, 2), (1, 3), (2, 3)}, Crit(i, X, Y , Q k p , Q k q ) end for Get (k * , p * , q * ) = arg max k∈[1:S K Y 3 ],(p,q)∈I,γ k ≥γ Crit(i, X, Y , Q k p , Q k q ) return Q k * 1 , Q k * 2 , Q k * 3 and Crit(i, X, Y , Q k * p * , Q k * q * )

Improved algorithm for first order indices based on single cluster membership

The previous algorithms are computationally intensive because of the exhaustive search over all partitions (2-partitions when optimizing indices of cluster membership, 3-partitions when optimizing indices of cluster membership differences). The associated costs are driven by the loop over these sets and is exponential with respect to the number of clusters of the pre-clustering set (O(2 K Y ), O(3 K Y )). In practice, it becomes intractable to consider K Y larger than 15 -20, which limits the exploration to large-sized clusters. This limit motivates the proposition of an improved algorithm that would take advantage of the properties of the cost functions (which were not taken into account in the different algorithms presented in the previous section). As the study then becomes dependent on the chosen criterion, we present here an improvement of the optimization based on SI C i (first order Sobol' index on single cluster membership).

We start by giving a property that helps building partitions from a set of elementary clusters (C k ) of the output space.

An aggregating property of clusters

We consider a model input X and suppose that it has a uniform distribution over an interval I. Let (C k ) k∈[1,K] be a set of non overlapping clusters of the output space. We consider the distribution of X given Y ∈ C k and denote as h C k the associated histograms obtained when discretizing I into n x bins. Let N be the total number of simulated points.

Then the discrete approximation of criterion SI C k , noted SI C k , can be written using h C k : This property states that elementary clusters with perfectly correlated histograms can not be in different sets of the optimal partition.

SI C k = n x nx j=1 h C k j (N - nx j=1 h C k j ) nx i=1 (h C k i - 1 n x nx j=1 h C k j ) 2

Algorithm

We derived an efficient algorithm named iRSA SM' by taking advantage of Proposition 1. The principle of iRSA SM' is to keep the overall structure of Algorithm iRSA SM and to improve the pre-processing step by using Proposition 1. More precisely, Proposition 1 is used to propose a merging of elementary clusters based on the correlation of their histograms. Indeed, if the histograms are perfectly correlated, then from Proposition 1, the associated elementary clusters can be merged.

The pre-clustering step of Algorithm iRSA SM' consists in two steps. Firstly, as for iRSA SM, a pre-clustering of Y using ClustF un Y is performed to get K Y clusters, but this time with a larger K Y (500-1000). Secondly, a clustering of these K Y clusters into K H meta-clusters is done using a distance-based method denoted ClustF un H . ClustF un H can be the same as ClustF un Y but is based on a different distance: ClustF un H works with an histogram-based correlation distance (e.g. d(h, h ) = 1 -Cor(h, h )) and is applied on the elementary histograms associated to the conditional distributions of X i given Y ∈ C k . The last step of the algorithm, the exhaustive search for the best partition over all possible merging of elementary clusters, is kept unchanged.

Compared to iRSA SM, Algorithm iRSA SM' has three additional parameters: the number n X of bins for the histograms computations, the generic method ClustF un H and the number of meta-clusters K H .

Algorithm iRSA SM'(i, X, Y , Crit, ClustF un Y , K Y , n x , ClustF un H , K H , γ) Apply ClustF un Y (K Y ) on Y to get K Y clusters (C 1 , .., C K Y ), with K Y large Compute all elementary histograms (h C k ) k=1..K Y for input X i discretized with n x bins Apply ClustF un H (K H ) on (h C1 , .., h C K Y ) to get K H meta-clusters ( Ĉ1 , .., ĈK H ) based on histogram correlation for all 2-partitions (P k , P k ) k∈[1:S K H 2 ] of ( Ĉ1 , .., ĈK H ) do Compute γ k = min Card(P k ), Card(P k ) Compute Crit(i, X, Y , P k ) end for Get k * = arg max k∈[1:S K H 2 ],γ k ≥γ Crit(i, X, Y , P k ) return (P k * , P k * ) and Crit(i, X, Y , P k * )
Note that this time, the algorithm is strongly dedicated to the optimization related to a single input X i due to the histogram computation step . It is thus only possible to re-use the first pre-clustering when performing the optimization for another input.

Application of the method on various models

Test model with 1d outputs: analytical resolution

In this section, we present the analytical resolution of the optimization problem on a simple 1D model in order to highlight the main issues of the iRSA approach.

Model and input distributions

We consider the model with two input factors X 1 and X 2 defined by :

Y = Sign(X 1 ) • |X 2 |
The model is studied over the domain defined by X 1 and X 2 having independent uniform distributions on [-1, 1]. Y thus also takes its values within [-1, 1].

Optimization problem

We consider the problem of finding for X 1 (resp. X 2 ) the partition (C, C) of the one-dimensional output space [-1, 1] the most influenced by the variations of X 1 (resp. X 2 ).

C * i = arg max C⊂[-1,1] SI C i , i = 1, 2
We limit the search to partitions that can be parameterized with a single cutting value y c or with two cutting values y c1 and y c2 . In the first case, it amounts to consider partitions of type 'A-B', and in the second case to partitions of type 'A-B-A' of the 1-dimensional output space. 7

SI C i (y c ) SI C i (y c1 , y c2 ) X 1 (1-|yc|) 2 1-y 2 c (|yc 2 |-|yc 1 |) 2 (yc 2 -yc 1 )(2-yc 2 +yc 1 ) X 2 |yc| 1+|yc| 1-yc 2 +yc1 2-yc 2 +yc 1 if y c2 ≤ 0 or y c1 ≥ 0 1-yc 2 +yc 1 2-yc 2 +yc 1 + 2 min(|yc 1 |,yc 2 ) (yc 2 -yc 1 )(2-yc 2 +yc 1 ) if y c1 ≤ 0 ≤ y c2

Analytical expressions of the region-based indices

The simplicity of this model allows to derive analytical expressions of the sensitivity indices (SI C i ) i=1,2 for the different types of partition using one cutting values y c or two cutting values y c1 and y c2 (see Appendix A for details on how to obtain these expressions).

Optimization

The expressions presented in Table 1 allow to derive the optimal cutting values that maximize the sensitivity indices for the two parameters in the different cases.

Single cutting value y c . In this case, as can be seen in Figure 1, it is easy to show that:

• SI C
1 is maximum for y c = 0. For this partition, we have

SI C 1 = 1 and SI C 2 = 0. • SI C
2 is maximum when y c reach the boundaries of the domain (y c = -1 or y c = 1). In this case, SI C 2 = 1 2 . For these partitions, we have SI C 1 = 0.

We conclude that the single cutting value case allows to find an optimal partition (C * 1 , C * 1 ) = ([-1, 0], [0, 1]) for input X 1 independently of the hypothesis on the number of cutting values (as SI C 1 maximum value is 1). However, the case of a single cutting value does not allow to find a partition fully explained by the variations of X 2 .

Two cutting values. In this case, as can be seen in Figure 2, we can show that:

• SI C
1 is maximum for (y c1 , y c2 ) = (-1, 0) and (y c1 , y c2 ) = (0, 1). This values are degenerated situations and correspond to the single cutting value case.

• SI C

2 is maximum when y c1 = -y c2 . For these partitions, we have SI C 2 = 1 and

SI C 2 = 0.
Considering the case of two cutting values allows to find a set of partitions with SI C 2 = 1. These partitions can be written

(C * 2 , C * 2 ) = ([-δ, δ], [-1, -δ] ∪ [δ, 1]), with δ ∈]0, 1[.

Conclusion of the example

Using the analytical expressions of the sensitivity indices given previously, the iRSA sensitivity approach allows to conclude that X 1 has the effect of leading negative output to positive outputs, and that X 2 has the effect of leading the output from the center of the domain toward its boundaries.

This example also highlights some issues about iRSA: • The optimal partition is associated to a model input and may be different depending on the input factor considered.

• Optimal partitions can be found even if the space of partitions is not completely explored. In this example, exploring partitions defined by one or two cutting values was sufficient to find global optima.

• We can have SI C * i > SI i (SI i being the standard first order Sobol' index of input factor X i ). In this example, we have SI C * 2 = 1 while SI 2 = 0: the information provided by iRSA is different from the one provided by a standard variance-based sensitivity analysis.

• The optimal partition may not be unique. In this example, we found a set of optimal partitions for X 2 .

Test model with 2D outputs: illustration of the different SI-based criteria 3.2.1. Model and input distributions

We consider the model having 2D outputs (Y 1 , Y 2 ) and four inputs (X 1 , X 2 , X 3 , X 4 ). Basically, the model draws a point in the plane at a distance related to X 4 and with an angle related to X 3 . This angle is defined from one of three given centers, with a choice of center related to X 1 and X 2 in interaction. More precisely, the model definition is the following:

Y 1 = c 1 + 0.4 • X 3 4 • cos(2πX 3 ) Y 2 = c 2 + 0.4 • X 3 4 • sin(2πX 3 ) where (c 1 , c 2 ) =     
(0.5, 0.25), if X 1 < 0.5 (0.25, 0.75), if X 1 ≥ 0.5 and X 2 < 0.5 (0.75, 0.75), if X 1 ≥ 0.5 and X 2 ≥ 0.5 All inputs have independent uniform distribution in [0, 1]. The sensitivity indices were computed using the SobolSalt function of the R package sensitivity [START_REF] Pujol | Sensitivity: Global sensitivity analysis of model outputs, R package[END_REF] with a sample size of n = 1000, producing 10000 output values.

Optimization setting

We illustrate the three numerical algorithms of Section 2 on this 2D toy model. More precisely, we apply for each input (X i ) i=1..4 Algorithm iRSA SM and its improved version iRSA SM' in order to find the clustering maximizing SI C i . We also apply Algorithm iRSA DM in order to find the clustering maximizing SI CC i . For algorithms iRSA SM and iRSA DM, we use a pre-clustering of the simulated results in K Y = 10 elementary clusters using a fuzzy k-means clustering approach (which is well adapted for small K Y ). This pre-clustering is presented in Figure 3. For Algorithm iRSA SM', the first pre-clustering was performed in K Y = 1000 elementary clusters. The histogram correlation method was then applied to obtain K H = 10 clusters before the application of the exhaustive search step with a size constraint parameter γ = 10%. A hierarchical clustering approach (more adapted for large K Y ) was used for ClustF un Y (with an euclidean distance on Y) and ClustF un H (with a correlation distance on histograms).

Method for displaying clusters

For Algorithm iRSA SM and iRSA DM, we present the clustering into two or three clusters and the associated region-based sensitivity indices. Algorithm iRSA SM' must be handled differently, as the algorithm does not guaranty clear boundaries between clusters. To handle this issue, we propose a different display of its results: we associate to each pixel a continuous color level (noted mmv, for mean membership value) linked to the proportion of output points belonging to each of the two clusters. Using this definition, pixels with output points from cluster A only (resp. B only) will be displayed in blue (resp. red), but pixels with output points from the two clusters will have a clearly different color. This display allows to quickly spot zones where the clusters boundaries are blurred and where interpreting the regions can be difficult. We present in Figure 4 the result obtained with Algorithm iRSA SM. For inputs X 1 and X 4 , the clustering found is associated to very high region-based indices (SI C * centers. Nevertheless, the optimized region-based sensitivity indices are very high. For input X 3 , the sensitivity score associated to the clustering is low and prevents any further interpretation. For X 2 , we obtain a relatively low sensitivity score (SI C * 2 = 0.27) but a clear structure of clusters: the optimal partition is made from points closest to top right center. This was expected, as X 2 drives the choice of the two centers of the top. But as this choice is made conditionally to X 1 , there is a strong interaction effect with X 1 that prevents from getting first order indices close to one.

Clustering based on SI CC with Algorithm iRSA DM

The optimization of criterion SI CC tries to find two regions whose transition from the one to the other is the most influenced by the variations of an input. By construction it leads to higher scores compared to optimizing SI C since SI C is a particular case of SI CC with C = ∅. When looking at Figure 5, we can verify this property for X 1 and X 4 . We can also see that the result obtained on X 3 still prevents any interpretation. For X 4 , the indices are still high, but the interpretation is not made easier. On the other hand for X 2 , the algorithm found an improved value of the first order index ( SI CC * 2 = 0.41 ) associated with a value of its total index TSI CC * = 0.97 that only X 2 reaches. This means that X 2 is the key model input that explains the transitions between the red cluster (points closest to the top right center) and the blue one (points closest to the top right center). Note again that it is not surprising to still get a strong interaction with X 1 given the model definition. This example shows the interest of considering the SI CC criterion. The results of the application of Algorithm iRSA SM' are presented in Figure 6. The algorithm allows to generate clusters with more flexibility in terms of boundaries and spatial resolution than Algorithm iRSA M. For X 1 and X 4 , the sensitivity indices were already well optimized using the first algorithm. However, this score is significantly improved for X 4 using Algorithm iRSA SM'. We can also note that the three red zones identified for X 4 have the clear expected circular shape, which was not possible to obtain with the first algorithm. For X 2 , the scores are still low (SI C * 3 = 0.30) but in this case, it is due to the definition of the model, not to a limit of the algorithm. A clear improvement was obtained for X 3 , for which the sensitivity score associated to the best clustering is strongly improved (SI C * 3 = 0.57 against SI C * 3 = 0.11 for Algorithm iRSA SM) and where the expected angular pattern starts to appear. The three red zones correspond approximately to the same angular section around their closest centers, which was expected given the model definition. To conclude, we have shown on this example that Algorithm iRSA SM' managed to automatically identify regions of the output space that maximize first-order region-based Sobol' indices. These regions were quite precisely delimited and are consistent with what was expected, given the model definition.

Test model producing time series: application of Algorithm iRSA SM' 3.3.1. Model and input distributions

We apply the iRSA approach on the environmental model CANTIS [START_REF] Garnier | Evaluation of a nitrogen transport and transformation model in a bare soil[END_REF], which simulates the decomposition of organic biomass in crop residues over time. We use the same model setting as in [START_REF] Roux | Cluster-based gsa: Global sensitivity analysis of models with temporal or spatial outputs using clustering[END_REF] and invite the reader to refer to this article to get more details. The model setting consists in ten uncertain model inputs and the analysis of one timevarying output: the zymogenous microbial biomass, here-after referred to as ZYB. The sensitivity indices are computed using the SobolSalt function of the R package sensitivity [START_REF] Pujol | Sensitivity: Global sensitivity analysis of model outputs, R package[END_REF] with a sample size n = 2500 producing 55000 output curves of ZYB. input

X 1 X 2 X 3 X 4 X 5 X 6 X 7 X 8 X 9 X 10 SI C * 1
0.056 0.054 0.057 0.057 0.042 0.049 0.038 0.843 0.868 0.031

Optimization setting

Algorithm iRSA M' is applied for each input (X i ) i=1..10 in order to find the clustering maximizing SI C i . A hierarchical clustering approach (well adapted to large K Y ) is used for ClustF un Y (with an euclidean distance on Y) and ClustF un H (with a correlation distance on histograms). The first pre-clustering is performed with K Y = 500 clusters and the second with K H = 10 clusters before the application of the exhaustive step with a size constraint of γ = 10%.

Method for displaying clusters

In this case, iRSA produces two sets of curves for each input factor. As in the 2D example, we represent the result using mean membership maps (see Section 3.2.3). We also plot for each cluster the quantiles of the distributions obtained at each time step. In order to quickly assess whether at a given time, the distribution of values of the two clusters differ a lot or not, we add an horizontal colorbar corresponding to the Kolmogorov-Smirnov test statistic (maximum value of the difference between cumulated distributions).

Results

We present in Table 2 the sensitivity indices SI C * 1 obtained when applying Algorithm iRSA SM' for each input of the CANTIS model in the considered uncertainty setting. Only X 8 and X 9 exhibits large values that deserve further analysis. For these two inputs, we present in Figure 7 the clusters using the visualization approach detailed in Section 3.3.3, which is necessary to handle large-sized data set (55000 clustered curves).

For input X 8 , we obtain SI C * 8 = 0.843 for the optimized clustering. Looking at Kolmogorov-Smirnov test statistic in Figure 7, we see that the two clusters differ the most at the beginning of the simulation. More precisely, when looking at the mean membership map, we see that the red cluster is composed mostly of curves that start to increase and the blue one of curves that start to decrease. The two clusters then overlap rapidly. The application of the method let us conclude that X 8 drives the occurrence of biomass dynamics that start by increasing (or equivalently that start by decreasing).

For input X 9 , a very high score SI C * 9 = 0.868 is obtained using Algorithm iRSA SM'. For this input, when looking Kolmogorov-Smirnov test statistic, it appears that what best distinguishes the two clusters is the level of ZYB values over the entire simulated time period, and particularly the final value. We can conclude that X 9 drives the occurrence of dynamics with highest (or lowest) biomass levels all along the considered time scale.

Conclusion

In this work, we introduced a new sensitivity analysis approach (named iRSA) for exploring models with various output dimensions by proposing a new perspective on the concept of Regional Sensitivity Analysis. Instead of a priori defining target model behaviors, we automatically detect the ones whose occurrences are best explained by the variations of the different model inputs. This is a new way of exploring simulation models. The method was illustrated by analytically solving a simple scalar case. Different algorithms have been proposed to compute numerical solutions in the case of complex models as well as adapted graphical representations to interpret the results. The relevance of the results was checked on a model with 2D outputs. iRSA was also applied on an environmental model producing time series, showing that interpretable information can be produced by the method. Future work on iRSA will be focused on i) the extension of the proposed sensitivity criteria in order to focus on groups of inputs and interactions, ii) attempting to derive efficient algorithms for other criteria (as we did for the first order indices), iii) the application of the method on models with other kind of complex outputs, such as spatial ones.

Proposition 1 .

 1 Let's consider two elementary clusters C k and C k with perfectly correlated histograms h C k and h C k = θh C k . Let's denote (C * , C * ) the partition of the output space built from the (C k ) that maximizes SI C (discretized first order index based on single cluster membership). Then C k and C k belong both to C * or to C * . Proof. See Appendix B.

Figure 1 : 2 CFigure 2 :

 122 Figure 1: First order indices SI C i as a function of a single cutting value yc for the two inputs of model Y = Sign(X 1 ) • |X 2 |

Figure 3 :

 3 Figure 3: Clustering of the 2D test model outputs into K Y = 10 clusters

Figure 4 :

 4 Figure 4: Results obtained with the SI C criterion and Algorithm iRSA SM on the 2D toy model. First line: region-based sensitivity indices of the best clustering for each input. Second line: visualization of the clusters for each input.

Figure 5 :

 5 Figure 5: Results obtained with the SI CC criterion and Algorithm iRSA DM on the 2D toy model. First line: region-based sensitivity indices of the best clustering for each input. Second line: visualization of the clusters (red and blue with neutral class in yellow) for each input.

3. 2 . 6 .

 26 Clustering based on SI C with Algorithm iRSA SM'

Figure 6 :

 6 Figure 6: Results obtained with the SI C criterion and Algorithm iRSA SM' on the 2D toy model. First line: region-based sensitivity indices of the best clustering for each input. Second line: visualization of the clusters for each input using mean membership maps.

Figure 7 :

 7 Figure 7: Clustering obtained with Algorithm iRSA SM' on the CANTIS model for input X 8 (first line) and X 9 (second line). First column: region-based sensitivity indices of the best clustering; second column: mean membership maps; third and fourth columns: quantiles (0.05,0.25,0.50,0.75,0.95) of the two clusters of curves at each time. The horizontal colorbar stresses significant differences between the two distributions Y (t)|(Y ∈ C) and Y (t)|(Y ∈ C). : yellow (resp. black) correspond to a high (resp. low) value of the Kolmogorov-Smirnov (ks) test statistics applied to the two distributions.

Table 1 :

 1 Analytical expressions of the region-based sensitivity indices for model Y = Sign(X 1 ) • |X 2 |

Table 2 :

 2 Optimized region-based sensitivity indices obtained for each input X i .

= 0.86, SI C * 4 = 0.80). These scores reveal a strong effect of these inputs in the occurrence of the obtained clusters. For X 1 , two large clusters are obtained with a boundary approximately equal to the line y

= 0.5. For X 4 , the clustering distinguishes a region made of three zones (each one close to one of the three centers) from a region made of points far from all centers. Note that the pre-clustering does not allow to get perfect circles around the

We derive for this model the expressions of SI C i when X 1 and X 2 have independent uniform distributions over [-1, 1] for two parameterizations of an interval C ⊂ [-1, 1] : one using a single cutting value y c and another using two cutting values y c1 and y c2 .

A.1 Parameterization with one cutting value y c

In this case, we study f (X 1 , X 2 ) = 1 sign(X1)|X2|≤yc

Variance.

Parameterization with two cutting values y c1 and y c2

In this case, we study

If y c1 and y c2 have the same sign:

If y c1 and y c2 have opposite sign:

We recall the expression of the discretized version of the first order sensitivity-based criterion used to performed the clustering based on single membership functions for an input X whose conditional distribution given Y ∈ C has an n x bins histogram noted h C :

Let us consider a cluster C 0 to which we either add a set of points with histogram h to form a cluster C 1 or two set of points with respective histograms h and h = θh forming a cluster C 2 . The criteria for these three sets are SI C0 , SI C1 , SI

C2

. We suppose that adding to C 0 the set of points with histogram h improves the clustering (C 0 , C 0 ), so that SI C0 ≤ SI

C1

. Then the principle of the proof is to show that there is a greater improvement in adding also the set with a perfectly correlated histogram h , ie that SI C1 < SI

C2

. Finally, supposing now then that there are two set C k and C k with perfectly correlated histograms, there are necessarily of the same side of the optimal partition made by grouping elementary clusters, otherwise grouping C k and C k would strictly improve the criterion, which is not possible.

So we assume now that SI

We denote:

Assuming