R: Computation of Sobol' indices

sobol_indices {sensobol}

R Documentation

Computation of Sobol' indices

Description

It allows to compute Sobol' indices up to the fourth-order using state-of-the-art estimators.

Usage

sobol_indices(
  matrices = c("A", "B", "AB"),
  Y,
  N,
  params,
  first = "saltelli",
  total = "jansen",
  order = "first",
  boot = FALSE,
  R = NULL,
  parallel = "no",
  ncpus = 1,
  conf = 0.95,
  type = "norm"
)

Arguments

`matrices`	Character vector with the required matrices. The default is `matrices = c("A", "B", "AB")`. See `sobol_matrices`.
`Y`	Numeric vector with the model output obtained from the matrix created with `sobol_matrices`. The numeric vector should not contain any NA or NaN values.
`N`	Positive integer, the initial sample size of the base sample matrix created with `sobol_matrices`.
`params`	Character vector with the name of the model inputs.
`first`	Estimator to compute first-order indices. Options are: `first = "saltelli"` (Saltelli et al. 2010). `first = "jansen"` (Jansen 1999). `first = "sobol"` (Sobol' 1993). `first = "azzini"` (Azzini et al. 2020).
`total`	Estimator to compute total-order indices. Options are: `total = "jansen"` (Jansen 1999). `total = "sobol"` (Sobol' 2001). `total = "homma"` (Homma and Saltelli 1996). `total = "janon"` (Janon et al. 2014). `total = "glen"` (Glen and Isaacs 2012). `total = "azzini"` (Azzini et al. 2020). `total = "saltelli"` (Saltelli et al. 2008).
`order`	Whether to compute "first", "second", "third" or fourth-order Sobol' indices. Default is `order = "first"`.
`boot`	Logical. If TRUE, the function bootstraps the Sobol' indices. If FALSE, it provides point estimates. Default is `boot = FALSE`.
`R`	Positive integer, number of bootstrap replicas. Default is NULL.
`parallel`	The type of parallel operation to be used (if any). If missing, the default is taken from the option "boot.parallel" (and if that is not set, "no"). For more information, check the `parallel` option in the `boot` function of the `boot` package.
`ncpus`	Positive integer: number of processes to be used in parallel operation: typically one would chose this to the number of available CPUs. Check the `ncpus` option in the `boot` function of the `boot` package.
`conf`	Confidence interval if `boot = TRUE`. Number between 0 and 1. Default is `conf = 0.95`.
`type`	Method to compute the confidence interval if `boot = TRUE`. Default is "norm". Check the `type` option in the `boot` function of the `boot` package.

Details

Any first and total-order estimator can be combined with the appropriate sampling design. Check Table 3 of the vignette for a summary of all possible combinations, and Tables 1 and 2 for a mathematical description of the estimators. If the analyst mismatches estimators and sampling designs, the function will generate an error and urge to redefine the sample matrices or the estimators.

For all estimators except Azzini et al. (2020)'s and Janon et al. (2014)'s, sobol_indices() calculates the sample mean as

\hat{f}_0=\frac{1}{2N} \sum_{v=1}^{N}(f(\mathbf{A})_v + f(\mathbf{B})_v)\,,

where N is the row dimension of the base sample matrix, and the unconditional sample variance as

\hat{V}(y) = \frac{1}{2N-1} \sum{v=1}^{N} ((f(\mathbf{A})_v - \hat{f})^2 + (f(\mathbf{B})_v - \hat{f})^2)\,,

where f(\mathbf{A})_v (f(\mathbf{B})_v) indicates the model output y obtained after running the model f in the v-th row of the \mathbf{A} (\mathbf{B}) matrix.

For the Azzini estimator,

\hat{V}(y) = \sum_{v=1}^{N} (f(\mathbf{A})_v - f(\mathbf{B})_v)^2 + (f(\mathbf{B}_A^{(i)})_v - f(\mathbf{A}_B^{(i)})_v) ^ 2

and for the Janon estimator,

\hat{V}(y)=\frac{1}{N} \sum_{v=1}^{N} \frac{f(\mathbf{A})_v^2 + f(\mathbf{A}_B^{(i)})_v^2}{2}-f_0^2

where f(\mathbf{A}_B^{(i)})_v (f(\mathbf{B}_A^{(i)})_v) is the model output obtained after running the model f in the v-th row of an \mathbf{A}_B^{(i)})_v (\mathbf{B}_A^{(i)})_v) matrix, where all columns come from \mathbf{A} (\mathbf{B}) except the i-th, which comes from \mathbf{B} (\mathbf{A}).

Value

A sensobol object.

References

Azzini I, Mara T, Rosati R (2020). “Monte Carlo estimators of first-and total-orders Sobol' indices.” arXiv. 2006.08232, https://arxiv.org/abs/2006.08232.

Glen G, Isaacs K (2012). “Estimating Sobol sensitivity indices using correlations.” Environmental Modelling and Software, 37, 157–166. doi:10.1016/j.envsoft.2012.03.014.

Homma T, Saltelli A (1996). “Importance measures in global sensitivity analysis of nonlinear models.” Reliability Engineering and System Safety, 52, 1–17. doi:10.1016/0951-8320(96)00002-6.

Janon A, Klein T, Lagnoux A, Nodet M, Prieur C (2014). “Asymptotic normality and efficiency of two Sobol index estimators.” ESAIM: Probability and Statistics, 18(3), 342–364. doi:10.1051/ps/2013040.

Jansen M (1999). “Analysis of variance designs for model output.” Computer Physics Communications, 117(1), 35–43. doi:10.1016/S0010-4655(98)00154-4.

Saltelli A, Annoni P, Azzini I, Campolongo F, Ratto M, Tarantola S (2010). “Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index.” Computer Physics Communications, 181(2), 259–270. doi:10.1016/j.cpc.2009.09.018.

Saltelli A, Ratto M, Andres T, Campolongo F, Cariboni J, Gatelli D, Saisana M, Tarantola S (2008). Global Sensitivity Analysis. The Primer. John Wiley and Sons, Ltd, Chichester, UK. doi:10.1002/9780470725184.

Sobol' IM (1993). “Sensitivity analysis for nonlinear mathematical models.” Mathematical Modeling and Computational Experiment, 1(4), 407–414.

Sobol' IM (2001). “Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates.” Mathematics and Computers in Simulation, 55(1-3), 271–280. doi:10.1016/S0378-4754(00)00270-6.

Examples

# Define settings
N <- 1000; params <- paste("X", 1:3, sep = ""); R <- 10

# Create sample matrix
mat <- sobol_matrices(N = N, params = params)

# Compute Ishigami function
Y <- ishigami_Fun(mat)

# Compute and bootstrap Sobol' indices
ind <- sobol_indices(Y = Y, N = N, params = params, boot = TRUE, R = R)

[Package sensobol version 1.1.5 Index]