# Robustness of the Sobol' indices to marginal distribution uncertainty

Global sensitivity analysis (GSA) quantifies the influence of uncertain variables in a mathematical model. The Sobol' indices, a commonly used tool in GSA, seek to do this by attributing to each variable its relative contribution to the variance of the model output. In order to compute Sobol' indices, the user must specify a probability distribution for the uncertain variables. This distribution is typically unknown and must be chosen using limited data and/or knowledge. The usefulness of the Sobol' indices depends on their robustness to this distributional uncertainty. This article presents a novel method which uses "optimal perturbations" of the marginal probability density functions to analyze the robustness of the Sobol' indices. The method is illustrated through synthetic examples and a model for contaminant transport.

## Authors

• 7 publications
• 5 publications
03/29/2018

### Robustness of the Sobol' indices to distributional uncertainty

Global sensitivity analysis (GSA) is used to quantify the influence of u...
03/31/2022

### Test comparison for Sobol Indices over nested sets of variables

Sensitivity indices are commonly used to quantify the relative influence...
08/07/2020

### An information geometry approach for robustness analysis in uncertainty quantification of computer codes

Robustness analysis is an emerging field in the domain of uncertainty qu...
01/15/2018

### An Efficient Method for Uncertainty Propagation in Robust Software Performance Estimation

Software engineers often have to estimate the performance of a software ...
02/06/2018

### Quantifying dependencies for sensitivity analysis with multivariate input sample data

We present a novel method for quantifying dependencies in multivariate d...
06/12/2019

### Global optimization using Sobol indices

We propose and assess a new global (derivative-free) optimization algori...
05/27/2020

### Assessing variable activity for Bayesian regression trees

Bayesian Additive Regression Trees (BART) are non-parametric models that...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

Uncertainties are ubiquitous in mathematical models of complex processes. Understanding how each source of input uncertainty contributes to the model output uncertainty is critical to the development and use of such models. Global sensitivity analysis (GSA) addresses this challenge by considering the inputs of the mathematical model as random variables (or implicitly assuming they are uniformly distributed) and computing a measure of their relative importance in determining the model output. A common tool for this purpose is the Sobol’ indices

[19, 20, 17]. In their classical framework, one assumes that the model inputs are independent random variables and is some function mapping to a scalar valued quantity of interest (QoI), a function of the model output. The Sobol’ indices apportion how much of the variance of is contributed by each input.

Computing the Sobol’ indices requires the user to specify a probability distribution for , draw samples from this distribution, and evaluate at these samples. However, the distribution of

is typically estimated from limited data and/or knowledge and hence is uncertain itself. This article studies the robustness of the Sobol’ indices to perturbations in the marginal distributions of

, i.e. the distributions of , . This study is motivated by problems where the only knowledge of is a uncertain nominal estimate, for instance a least squares estimate with noisy data, and possibly loose upper and lower bounds. This is common in many applications. For such problems, the user often assumes that each follows a symmetric distribution centered at the nominal estimate with some variance, usually chosen with limited knowledge. An example is centering a uniform distribution at the nominal estimate. The approach proposed in this article may be applied by (i) giving each a distribution with support on a wide interval with greater probability near the nominal estimate, (ii) computing the Sobol’ indices using these distributions, (iii) and assessing the Sobol’ indices robustness with respect to changes in the distributions.

In the authors previous work [10], the robustness of Sobol’ indices was studied by taking perturbations of the joint probability density function (PDF) of . This approach was useful for studying robustness with respect to statistical dependencies in . This article focuses on cases where the inputs are independent and the uncertainty is exclusively in their marginal distributions. We modify the framework introduced in [10] to develop a novel method to study the robustness of Sobol’ indices with respect to changes in the marginal distributions. It is not a generalization, specialization, nor improvement of [10], but rather a novel method which leverages ideas from [10].

The robustness of Sobol’ indices with respect to changes in the marginal distributions of has been considered in the life cycle analysis literature [14] and ecology literature [16]. “Robust sensitivity indicators” are defined in [6] to account for uncertainty in the marginal distributions of . A generalization of the analysis of variance decomposition is considered in [2] to facilitate studying the dependence of Sobol’ indices on the distribution of . Other work such as [7, 8, 4] explore the robustness of Sobol’ indices to distributional uncertainty and [15, 11, 3, 1] explore more general questions of robustness for other computed quantities. In all of these cases, the user is required to either specify the parametric forms of distributions a priori, compute additional evaluations of beyond those needed to estimate Sobol’ indices, or both. The approach proposed in this article requires neither.

Section 2 reviews the Sobol’ indices and sets the notation for the article. We present the mathematical foundations of our proposed marginal PDF perturbation method in Section 3 and an algorithmic summary is given in Section 4. User guidelines are considered in Section 5. Numerical results are presented in Section 6. Concluding remarks are made in Section 7.

## 2 Review of Sobol’ indices

Let be a random variable supported on , , , and . Assume that the inputs are independent and is square integrable. For , let denote the corresponding subset of input variables. Then admits the following ANOVA (analysis of variance) decomposition,

 f(X)=f0+p∑i=1fi(Xi)+p∑i=1p∑j=i+1fi,j(Xi,Xj)+⋯+f1,2,…,p(X) (1)

where

 f0=E[f(X)], fi(Xi)=E[f(X)|Xi]−f0, fi,j(Xi,Xj)=E[f(X)|Xi,Xj]−fi(Xi)−fj(Xj)−f0, ⋮ fu(Xu)=E[f(X)|Xu]−∑v⊂ufv(Xv).

The independence of the inputs implies that

 E[fu(Xu)fv(Xv)]=0% whenu≠v. (2)

Computing the variance of both sides of (1) and applying (2) yields

 Var(f(X))=p∑i=1Var(fi(Xi))+p∑i=1p∑j=i+1Var(fi,j(Xi,Xj))+⋯+Var(f1,2,…,p(X)). (3)

Hence the variance of may be decomposed into contributions from each subset of input variables. The Sobol’ index for the variables is defined as

 Su=Var(fu(X))Var(f(X)),

and may be interpreted as the relative contribution of to the variance of . In practice, one frequently computes the first order Sobol’ indices , and the total Sobol’ indices defined by

 Tk=∑k∈uSu,k=1,2,…,p.

It easily observed that . The first order Sobol’ index may be interpreted as the relative contribution of variable to the variance of by itself and the total Sobol’ index may be interpreted as the relative contribution of to the variance of by itself and through interactions with other variables. The total Sobol’ index may also be defined for a subset of variables in an analogous fashion. This article will focus on the indices ; though the results are easily extended to general indices or , .

There are a plurality of ways to estimate the first order and total Sobol’ indices [17]. We focus on the Monte Carlo integration scheme presented in [17]. The indices and may be written probabilistically as

The properties of conditional expectation and Fubini’s theorem yield

 Var(E[f(X)|Xk])=E[f(X)(f(X′k,X∼k)−f(X′))],

and

 E[Var(f(X)|X∼k)]=12E[(f(X)−f(X′k,X∼k))2],

where is an independent copy of .

The indices may be estimated by Monte Carlo integration using evaluations of , where is the number of Monte Carlo samples. Specifically, the samples are generated by constructing 2 matrices of size , call them and , where the rows correspond to independent samples of . Then additional matrices of size , denote them as , , are generated by setting the column of equal to the column of , , and the column of equal to the column of . The set of indices may be estimated using evaluations of at the rows of , , and , . This scheme is presented in Lines 1-6 of Algorithm 1.

The reader is directed to [19, 20, 12, 9, 17, 21, 18, 13, 22], and references therein, for additional discussion about the Sobol’ indices.

## 3 Marginal PDF Perturbations

We study the robustness of the Sobol’ indices to distributional uncertainty by computing a derivative of the Sobol’ index with respect to the PDFs of . To this end we make the following assumptions:

1. is a compact interval, ,

2. admits a PDF , ,

3. is continuous on , ,

4. , ,

5. is bounded on .

Without loss of generality, under the assumptions above, assume that . Distributions with unbounded support may be truncated to satisfy these assumptions, hence our method is applicable on a large class of problems.

In order to perturb the PDFs while preserving the non-negativity and unit mass properties of PDFs, we define the Banach spaces , , as all bounded functions on with norm

 ∣∣∣∣ψi∣∣∣∣Vi=∣∣∣∣ψiϕi∣∣∣∣L∞(Ωi),

where the norm is the traditional supremum norm. The product Banach space is equipped with norm

 ∣∣∣∣(ψ1,ψ2,…,ψp)∣∣∣∣V=max1≤i≤p∣∣∣∣ψi∣∣∣∣Vi.

Then for with , , the functions

 ϕi+ψi1+∫Ωiψi(xi)dxi

are PDFs supported on , , i.e. are non-negative and integrate to 1 over .

For , define the first order and total Sobol’ indices as the operators , ,

 Sk(η)=Fk(η)Hk(η)andTk(η)=Gk(η)Hk(η),

where are the operators

 Fk(η)=1p∏i=1∫Ωiηi(y)dy∫Ω×Ωf(x)(f(x′k,x∼k)−f(x′))p∏i=1ηi(xi)ηi(x′i)dxdx′,
 Gk(η)=121∫Ωkηk(y)dy∫Ω×Ωk(f(x)−f(x′))2p∏i=1ηi(xi)ηk(x′k)dxdx′k,
 Hk(η)=∫Ωf(x)2p∏i=1ηi(xi)dx−1p∏i=1∫Ωiηi(y)dy(∫Ωf(x)p∏i=1ηi(xi)dx)2.

For with , and , , are the first order and total Sobol’ indices computed when has the joint PDF

 p∏i=1ηip∏i=1∫Ωiηi(y)dy.

Theorem 3.1 gives the Fréchet derivative of the Sobol’ index with respect to the PDFs.

###### Theorem 3.1

The operators and , , are Fréchet differentiable at the nominal PDF and the Fréchet derivatives and , , are the bounded linear operators given by

 DSk(ϕ)ψ=DFk(ϕ)ψHk(ϕ)−Sk(ϕ)DHk(ϕ)ψHk(ϕ),DTk(ϕ)ψ=DGk(ϕ)ψHk(ϕ)−Tk(ϕ)DHk(ϕ)ψHk(ϕ),

where

 DFk(ϕ)ψ=∫Ω×Ωf(x)(f(x′k,x∼k)−f(x′))(p∑i=1ψi(xi)ϕi(xi)+ψi(x′i)ϕi(x′i))p∏i=1ϕi(xi)ϕi(x′i)dxdx′ −(∫Ω(p∑i=1ψi(xi)ϕi(xi))p∏i=1ϕi(xi)dx)∫Ω×Ωf(x)(f(x′k,x∼k)−f(x′))p∏i=1ϕi(xi)ϕi(x′i)dxdx′,
 DGk(ϕ)ψ= 12∫Ω×Ωk(f(x)−f(x′))2(ψk(x′k)ϕk(x′k)+p∑i=1ψi(xi)ϕi(xi))p∏i=1ϕi(xi)ϕk(x′k)dxdx′k − 12∫Ωkψk(y)dy∫Ω×Ωk(f(x)−f(x′))2p∏i=1ϕi(xi)ϕk(x′k)dxdx′k,
 DHk(ϕ)ψ= ∫Ωf(x)2(p∑i=1ψi(xi)ϕi(xi))p∏i=1ϕi(xi)dx +(∫Ω(p∑i=1ψi(xi)ϕi(xi))p∏i=1ϕi(xi)dx)(∫Ωf(x)p∏i=1ϕi(xi)dx)2 −2(∫Ωf(x)p∏i=1ϕi(xi)dx)(∫Ωf(x)(p∑i=1ψi(xi)ϕi(xi))p∏i=1ϕi(xi)dx).

###### Proof

One may easily observe that in a neighborhood of (assuming is not constant). It is sufficient to compute the Fréchet derivatives of , and as the Fréchet derivatives of and will follow from the quotient rule for differentiation. The Fréchet derivative of the operator

 (ϕ1,ϕ2,…,ϕp)↦∫Ωp∏i=1ϕi(xi)dx,

mapping from to , acting on , is easily shown to be

 ∫Ωp∑j=1ψj(xj)ϕj(xj)p∏i=1ϕi(xi)dx.

Similar results hold for the Fréchet derivatives of

 ∫Ωf(x)p∏i=1ϕi(xi)dxand∫Ωf(x)2p∏i=1ϕi(xi)dx.

The Fréchet derivative of

follows from the sum/difference, product, and chain rules for differentiation. Similar arguments may be used to compute Fréchet derivatives of

and .

Theorem 3.1 can be generalized to Sobol’ indices or for any ; we omit this generalization for simplicity. If the Sobol’ indices are estimated with the Monte Carlo integration scheme presented in Section 2, then the action of the Fréchet derivatives on any may be estimated using the existing evaluations of . Hence the additional computation needed to estimate the Fréchet derivative is negligible.

We would like determine an optimal perturbation of the marginal PDFs which yields the greatest change in the Sobol’ indices. For notational simplicity, let denote a generic Sobol’ index, which may be or for any . The locally optimal perturbation direction is given by the , , which maximizes the absolute value of the Fréchet derivative of . To estimate it, we discretize by partitioning each interval into disjoint subintervals, , i.e. , , and define basis functions as

 ψji(X)=χRji(Xi)ei,j=1,2,…,Mi,i=1,2,…,p

where is the indicator function of a set and is the

canonical unit vector in

. Then

 VM=span{ψji|j=1,2,…,Mi,i=1,2,…,p}⊂V,

is a subspace of on which we will perform our robustness analysis. Discretizing with piecewise constant functions may be viewed as approximating all bounded functions by discretizing them on a coarse grid which is constrained by the existing samples and evaluations of .

Then the locally optimal perturbation is solution of the optimization problem

 maxψ∈VM||ψ||V=1|DI(ϕ)ψ|. (4)

Exploiting properties of the basis we chose for , (4) may be solved in closed form yielding the optimal perturbation direction

 ψ⋆=(ψ⋆1,ψ⋆2,…,ψ⋆p)=p∑i=1M∑j=1ajiψji

where

 aji=(infxi∈Rjϕi(xi))%sign(DI(ϕ)ψji).

Because of the norm on , this formulation gives equal weight to each marginal distribution. In actuality, we would like to focus our effort on taking perturbations on the “most important” marginals. To this end, we define the weights

and define the perturbed PDF as

 p∏i=1ϕi(xi)+δδiψ⋆i1+δδi∫Ωiψ⋆i(xi)dxi, (5)

where , ; the determination of will be discussed in Section 4. The perturbed PDF depends upon the partitions , , and perturbation size , which depends on user specified parameters, see Section 4.

We will refer to the distribution of as the nominal distribution and the distribution with PDF (5) as the perturbed distribution. Likewise, we refer to the Sobol’ indices computed using the nominal distribution as nominal Sobol’ indices. We compute perturbed Sobol’ indices, where the inputs follow the perturbed distribution, by reweighing the evaluations of used to estimate the nominal Sobol’ indices. Hence the computational effort to compute perturbed Sobol’ indices is negligible.

## 4 Algorithm

Evaluating the perturbed PDF requires estimation of the Sobol’ index Fréchet derivative to compute , and the determination of . As previously discussed, the Sobol’ index Fréchet derivative is easily estimated using existing evaluations of . The scalar may be determined by finding the largest and smallest possible values in

for which the perturbed index estimator’s sample standard deviation is sufficiently small. Specifically, a threshold should be defined by the user and a value of

is accepted if the ratio of sample standard deviations for the perturbed indices and nominal indices is below the threshold.

Algorithm 1 provides an overview of the Sobol’ index computation with marginal distribution robustness post processing. The algorithm requires four inputs:

1. , the number of Monte Carlo samples

2. , a partition of ,

3. , the threshold for admissible increases in the Sobol’ index estimator’s sample standard deviation

4. , an integer used to compute admissible value of

The determination of is problem dependent and must be done by the user to compute Sobol’ indices.

The partitions , , may be a chosen in many ways. A general partitioning strategy is to specific a number of subintervals and partition each into

subintervals according to the quantiles of the marginal distribution of

. Alternatively, the user may use problem specific information to customize the partition, for instance, taking a finer partition in regions of greater uncertainty or adapting based on the uncertainty in each marginal.

The threshold bounds the increase in estimation error introduced by reweighing samples to compute perturbed Sobol’ indices. If the original Sobol’ index estimation has small estimator standard deviation, then is considered a reasonable threshold to permit nontrivial perturbations while preserving reasonable accuracy of the estimator.

To determine admissible values of , we consider a function which inputs a perturbation size and returns the maximum (over all indices) ratio of sample standard deviations of the perturbed and nominal Sobol’ index estimators. We evaluate by subsampling the data and computing a collection of estimators and their standard deviations. Finding admissible values of is equivalent to finding solutions of . Intuition (and direct calculations) indicates that is approximately a quadratic function centered at . Admissible values of are determined by evaluating at equally spaced points in and accepting those for which . Since it is discretizing an interval, does not need to be very large, is suggested as a reasonable choice to balance the computational cost (the for loop in Line 14 of Algorithm 1) and the necessary discretization to determine admissible values of . The evaluations of at points may be plotted to confirm that is sufficient.

Our choices for the inputs described above are based on numerical experimentation, intuition, and experience from the authors previous work in [10].

Algorithm 1 returns the nominal Sobol’ indices and a collection of perturbed Sobol’ indices. The user may query these perturbed Sobol’ indices and visualize measures of robustness alongside the indices themselves. In addition, samples may be drawn from the perturbed distributions. The ratios of the nominal PDF and perturbed PDFs are computed (though not explicitly written) in Algorithm 1. Using candidate samples from the nominal distribution and the ratio of PDFs provides an easy way to generate samples from the perturbed distributions with rejection sampling. Visualizing the perturbed marginal distributions provides additional insights into the regions of parameter space which cause a lack of robustness. The examples in Section 6 illustrate efficient ways to visualize the output of Algorithm 1 and make inferences with it.

## 5 User guidelines

The method presented in Section 3 only considers perturbations of the one dimensional marginals PDFs, we refer to this as the marginal perturbation approach. In the authors previous work [10], perturbation were taken on the joint PDF, we refer to this as the joint perturbation approach. Other forms of perturbations may be considered as well. For instance, assuming a block dependency structure motivates a form of PDF perturbations on the marginal PDFs for each block. This approach would be a hybrid between the marginal perturbation approach of this article and joint perturbation approach of [10]. In this section, we compare and contrast important features of the marginal perturbation and joint perturbation approaches to provide some principles for a user to consider when deciding which method is appropriate for a given application. Since the two approaches seek to solve different problems, the user must be cognizant of their differences when determining which to use in an application.

If the inputs are assumed to be independent for lack of better knowledge, but dependencies are suspected, then the joint perturbation approach is more appropriate. If the inputs are believed to be independent, the joint perturbation approach may give pessimistic results by introducing dependencies which are not realistic, hence the marginal perturbation approach is superior.

In the joint perturbation approach, the user must generate a partition of which ensures a sufficient number of samples in each subset. If (the number of uncertain inputs) is large then the joint perturbation approach yields a partition which is coarse in each input domain . This problem is amplified if the nominal distribution of is far from being uniform. In contrast, the partitioning is trivial in the marginal perturbation approach because it only requires partitioning each , , separately. Hence it allows a richer class of perturbations in the marginal distributions. In addition, the user may insert a priori knowledge about each marginal distribution separately, for instance, in the partitions input to Algorithm 1.

In the marginal perturbation approach, each marginal is weighted by , a measure of the Fréchet derivative’s magnitude. This weighting strategy manifests itself differently in the joint perturbation approach, where the input variables are prioritized by the partitioning algorithm. The optimality of the weighting scheme is unclear in both cases; however, the simplicity of the weighting scheme for the marginal perturbation approach makes it more appealing.

The user must prioritize whether there is greater uncertainty in the marginal distributions or the dependency structure. The number of uncertain inputs and characteristics of their nominal distribution influence the performance of the two approaches and must also be considered. Subsection 6.1 gives a numerical illustration comparing the two approaches.

## 6 Numerical results

Three examples are given below. The first example compares the marginal perturbation approach of this article with the joint perturbation approach of [10]. The second is a synthetic example serving to illustrate our marginal perturbation method. The third is an application of robustness analysis to an advection diffusion contaminant transport problem. For simplicity and conciseness, we only report results for the total Sobol’ indices.

### 6.1 Comparison of marginal perturbation and joint perturbation approaches

Let

 f(X)=10∑i=1aiXi, (6)

where , . We consider two cases, first when has a uniform distribution on the interval , , and second, when

has a truncated normal distribution on the interval

, . In the truncated normal case, the mean of is and its standard deviation is , .

The marginal perturbation and joint perturbation approaches are compared by applying them to analyze the robustness of the total Sobol’ indices of (6), computed using Monte Carlo samples, for each of the two cases. Figure 1 displays the nominal, maximum, and minimum Total Sobol’ indices, where the maximum and minimum are determined by searching all perturbed total Sobol’ indices with . The top of the green bar and bottom of the magenta bar denote the maximum and minimum total Sobol’ index, respectively, as determined by the marginal perturbation approach; the blue region between them denotes the nominal total Sobol’ index. The top of the orange bar and bottom of the yellow bar denote the maximum and minimum total Sobol’ index, respectively, as determined by the joint perturbation approach; likewise the blue region between them denotes the nominal total Sobol’ index. The left panel corresponds to the uniform marginals case and the right panel corresponds to the truncated normal marginals case.

In both cases, the marginal perturbation approach found larger perturbation of the total Sobol’ indices than the joint perturbation approach, i.e. the marginal perturbation approach discovered a greater lack of robustness than the joint perturbation approach. This is because it was able to take a finer partition of and hence find “better” perturbations.

The total Sobol’ indices in truncated normal marginals case are more robust than in the uniform marginals case. This occurs because the largest perturbations of corresponds to putting greater probability near the boundaries, i.e. for near 0 and 1. The truncated normal marginals case puts less probability near the boundaries so perturbations of this nominal PDF will not change the total Sobol’ indices as much as in the uniform marginals case. Also notice that the joint perturbation approach yields very small perturbation for the truncated marginals case. This is because there are low probability tails in several marginal distributions which hampers the partitioning algorithm. The marginal perturbation approach does not suffer this limitation.

Using the marginal perturbation approach, the maximum and minimum perturbed total Sobol’ indices are approximately symmetric about the nominal total Sobol’ indices. In the joint perturbation approach, statistical dependencies are introduced which result in a decrease of the magnitude of the total Sobol’ indices, see [9, 10]

. As a result, the joint perturbation approach yields a skewness where the largest perturbations correspond to decreasing the total Sobol’ indices.

### 6.2 Synthetic example

To emulate a problem where there is uncertainty in the support of , let be uniformly distributed on the interval , where and are themselves random variables which are uniformly distributed on and , respectively, . Consider defined by

 f(X)=2X2e−2X1+X23. (7)

It is apparent from the structure of (7) that will vary more in the direction when is close to 0 and less when is close to 1. We study the sensitivity of to , and its robustness with respect to changes in the marginal distributions. Figure 2 displays the total Sobol’ indices of (7) in blue, along with the total Sobol’ indices under perturbations of the marginal distributions which maximize in red and minimize in cyan.

Figure 3 displays samples from the nominal distribution in blue (left), the perturbed distribution which maximizes in red (center), and the perturbed distribution which minimizes in cyan (right). In each column, the top row is the marginal distribution for , the middle row is the marginal distribution for , and the bottom row is the marginal distribution for . Notice that, as expected, the distribution of is perturbed to give greater probability near 0 (1) to maximize (minimize) . Further, to maximize (minimize) , the distribution of is concentrated toward the boundary (interior) which has the effect of maximizing (minimizing) the variance of . Conversely, the distribution of is concentrated toward the interior (boundary) to minimize (maximize) the variance of , since decreasing (increasing) the variance of will increase (decrease) .

In practice, the support of is frequently unknown and in many cases chosen ad hoc. Using a mixture distribution to model and coupling it with marginal perturbation robustness analysis, is a simple way to account for this uncertainty. Further, visualizing the perturbed distributions indicates whether the support will significantly influence the Sobol’ indices.

In this example, we analyze the sensitivity of a contaminant concentration to uncertainties in its source and dynamics. Specifically, let be the steady state contaminant concentration at a given point in the domain and our QoI be the average concentration in the subregion . To evaluate the QoI , we solve

 −ϵ△u+v⋅∇u=s on (0,1)2 u=0 on {0}×(0,1) u=0 on (0,1)×{0} n⋅∇u=ν1u on {1}×(0,1) n⋅∇u=ν2u on (0,1)×{1}

where the velocity field and source are given by

 v(y1,y2)=(α1(y1+0.5),α2(y2+0.5))ands(y1,y2)=βe−γ((y1−ξ1)2+(y2−ξ2)2),

and

denotes the outer normal vector to the boundary. The partial differential equation is solved using Chebfun

[5]. The QoI is then computed by integrating over . In this example,

 X=(ϵ,α1,α2,ξ1,ξ2,γ,β,ν1,ν2)

is a vector of 9 uncertain inputs which are assumed to be independent uniform random variables with uncertainty about their nominal values, given in Table 1.

Figure 4 illustrates several features of this advection diffusion example. The left panel of Figure 4 shows contours of the average source (computed from 10,000 samples), with the average velocity field overlaid with red arrows. The right panel of Figure 4 displays contours of the average solution . In both panels, cyan lines indicate the Robin boundaries and magenta lines indicate the Dirichlet boundaries and a black box encloses the region over which we integrate to evaluate .

The total Sobol’ indices of are computed via Monte Carlo integration with samples, followed by marginal distribution robustness analysis post processing. In a similar fashion to Figure 1, Figure 5 displays the nominal, maximum, and minimum total Sobol’ indices observed over all perturbations with . The top of the green bar and bottom of the magenta bar denote the maximum and minimum values, respectively; the blue region between them denotes the value of the nominal total Sobol’ index. The green and magenta bars in Figure 5 do not correspond to a particular set of Sobol’ indices for a fixed perturbation, but rather the extreme values observed over all perturbations. We observe that , , and exhibit minimal influence on our QoI and do not change significantly when the distribution of is perturbed. The source magnitude is found to be the most important parameter (under the nominal distribution), but its total Sobol’ index varies significantly when the marginal distributions are perturbed.

To further explore the influence of , we analyze the particular perturbations which cause the greatest change in its total Sobol’ index. Figure 6 displays the nominal total Sobol’ indices in blue, the perturbed total Sobol’ indices which maximize the total Sobol’ index for in red, and the perturbed total Sobol’ indices which minimize the total Sobol’ index for in cyan. This is complemented by Figure 7 which shows the marginal distributions for each parameter using the nominal distribution in blue, the perturbed distribution which maximize the total Sobol’ index for in red, and the perturbed distribution which minimize the total Sobol’ index for in cyan. We observe that in the maximum case (red), the total Sobol’ index for is much larger than all other total Sobol’ indices. The perturbed distribution in this case exhibit several features:

1. greater probability near the boundaries of ’s support,

2. greater probability near the interior of ’s support,

3. greater probability for smaller values of and ,

4. greater probability for larger values of and .

The minimum case (cyan) gives a different conclusion where is the fifth most important parameter. The features of the perturbed distribution in this case are opposite of the maximum case (red). In particular, the source location and advection velocity become more important by giving greater probability to the source being in the top left region of the domain (smaller , larger ), and having having greater probability in the tails of the horizontal advection speed .

## 7 Conclusion

A novel method has been presented which considers perturbations of the one dimensional marginal PDFs to measure the robustness of Sobol’ indices to distributional uncertainty. In its construction, this method assumes that the inputs are independent and the uncertainty is strictly in their one dimensional marginal distributions. This is a common occurrence in practice when, for instance, application specific information indicates that the inputs do not posses statistical dependences, but their statistical characterization is poorly known because of sparse and/or noisy data. In contrast, the authors previous work [10] provides a method to study robustness with respect to changes in the joint PDF, thus accounting for possible statistical dependencies. The source of distributional uncertainty is problem dependent and must be determined by the user. Since both methods are post processing steps after computing Sobol’ indices, the user may apply both methods and compare their results. Future work may include exportation of other PDF perturbation forms such as block dependency structures, and extensions of the PDF perturbation approach to other global sensitivity analysis methods.

## Software Availability

Matlab codes are available at

which implement the joint perturbation approach of [10] and the marginal perturbation approach presented in this article.

## References

• [1] Richard J. Beckman and Michael D. McKay. Monte Carlo estimation under different distributions using the same simulation. Technometrics, 29(2):153–160, 1987.
• [2] Emanuele Borgonovo, Max D. Morris, and Elmar Plischke. Functional ANOVA with multiple distributions: Implications for the sensitivity analysis of computer experiments. SIAM/ASA J. Uncertain. Quantif., 6(1):397–427, 2018.
• [3] Stephen E. Chick. Input distribution selection for simulation experiments: Accounting for input uncertainty. Operations Research, 49(5):744–758, 2001.
• [4] Areski Cousin, Alexandre Janon, Veronique Maume Deschamps, and Ibrahima Niang. On the consistency of Sobol’ indices with respect to stochastic ordering of model parameters.
• [5] Tobin A. Driscoll, Nicholas Hale, and Lloyd N. Trefethen, editors. Chebfun Guide. Pafnuty Publications, Oxford, 2014.
• [6] Lei Gao, Brett A. Bryan, Martin Nolan, Jeffery D. Connor, Xiaodong Song, and Gang Zhao. Robust global sensitivity analysis under deep uncertainty via scenario analysis. Environmental Modelling & Software Software, 76:154–166, 2016.
• [7] Evelyne A Groena and Reinout Heijungs. Ignoring correlation in uncertainty and sensitivity analysis in life cycle assessment: what is the risk? Environmental Impact Assessment Review, 62:98–109, 2017.
• [8] Jim Hall. Uncertainty-based sensitivity indices for imprecise probability distributions. Reliability Engineering and System Safety, 91:1443–1451, 2006.
• [9] Joseph Hart and Pierre Gremaud. An approximation theoretic perspective of Sobol’ indices with dependent variables. Internation Journal for Uncertainty Quantification, 8(6):483–493, 2018.
• [10] Joseph Hart and Pierre Gremaud. Robustness of the sobol’ indices to distributional uncertainty. https://arxiv.org/abs/1803.11249. Under review, 2018.
• [11] Zhaolin Hu, Jing Cao, and L. Jeff Hong. Robust simulation of global warming policies using the dice model. Management Science, 58(12):2190–2206, 2012.
• [12] Bertrand Iooss and Paul Lemaître. A review on global analysis methods. In G. Dellino and C. Meloni, editors, Uncertainty management in simulation-optimization of complex systems, chapter 5, pages 543–501. Springer, 2015.
• [13] Sergei Kucherenko, Stefano Tarantola, and Paola Annoni. Estimation of global sensitivity indices for models with dependent variables. Computer Physics Communications, 183:937–946, 2012.
• [14] Martino Lacirignola, Philippe Blanc, Robin Girard, Paula Pérez-López, and Isabelle Blanc. LCA of emerging technologies: addressing high uncertainty on inputs’ variability when performing global sensitivity analysis. Science of the Total Environment, 578:268–280, 2017.
• [15] Antony Millner, Simon Dietz, and Geoffrey Heal. Scientific ambiguity and climate policy. Environ Resource Econ, 55:21–46, 2013.
• [16] Livia Paleari and Roberto Confalonieri. Sensitivity analysis of a sensitivity analysis: We are likely overlooking the impact of distributional assumptions. Ecological Modelling, 340:57–63, 2016.
• [17] Clémentine Prieur and Stefano Tarantola. Variance-based sensitivity analysis: Theory and estimation algorithms. In Roger Ghanem, David Higdon, and Houman Owhadi, editors, Handbook of Uncertainty Quantification. Springer, 2017.
• [18] Andrea Saltelli, Paola Annoni, Ivano Azzini, Francesca Campolongo, Marco Ratto, and Stefano Tarantola. Variance based sensitivity analysis of model output. design and estimator for the total sensitivity index. Computer Physics Communications, 181:259–270, 2010.
• [19] Ilya M. Sobol’. Sensitivity estimates for non linear mathematical models. Math. Mod. Comp. Exp., 1:407–414, 1993.
• [20] Ilya M. Sobol’. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Mathematics and Computers in Simulation, 55:271–280, 2001.
• [21] Ilya M. Sobol’. Theorems and examples on high dimensional model representation. Reliability Eng. Sys. Safety, 79:187–193, 2003.
• [22] Stefano Tarantola and Thierry A. Mara. Variance-based sensitivity indices of computer models with dependent inputs: The Fourier Amplitude Sensitivity Test. International Journal for Uncertainty Quantication, 7(6):511–523, 2017.