sensitivity (version 1.14.0)

sobolroalhs: Sobol' Indices Estimation Using Replicated OA-based LHS

Description

sobolroalhs implements the estimation of the Sobol' sensitivity indices introduced by Tissot & Prieur (2015) using two replicated designs (Latin hypercubes or orthogonal arrays). This function estimates either all first-order indices or all closed second-order indices at a total cost of $2*N$ model evaluations. For closed second-order indices $N=q^2$ where $q >= d-1$ is a prime number corresponding to the number of levels of the orthogonal array, and where $d$ indicates the number of factors.

Usage

sobolroalhs(model = NULL, factors, N, p=1, order, tail=TRUE, conf=0.95, nboot=0, ...)
"tell"(x, y = NULL, ...)
"print"(x, ...)
"plot"(x, ylim = c(0,1), ...)

Arguments

model

a function, or a model with a predict method, defining the model to analyze.

factors

an integer giving the number of factors, or a vector of character strings giving their names.

N

an integer giving the size of each replicated design (for a total of $2*N$ model evaluations).

p

an integer giving the number of model outputs.

order

an integer giving the order of the indices (1 or 2).

tail

a boolean specifying the method used to choose the number of levels of the orthogonal array (see "Warning messages").

conf

the confidence level for confidence intervals.

nboot

the number of bootstrap replicates.

x

a list of class "sobolroalhs" storing the state of the sensitivity study (parameters, data, estimates).

y

a vector of model responses.

ylim

y-coordinate plotting limits.

...

any other arguments for model which are passed unchanged each time it is called.

Value

sobolroalhs returns a list of class "sobolroalhs", containing all the input arguments detailed before, plus the following components:

Warning messages

Details

sobolroalhs automatically assigns a uniform distribution on [0,1] to each input. Transformations of distributions (between U[0,1] and the wanted distribution) have to be realized before the call to tell() (see "Examples").

Missing values (i.e NA values) in outputs are automatically handled by the function.

This function also supports multidimensional outputs (matrices in y or as output of model). In this case, aggregated Sobol' indices are returned (see sobolMultOut).

References

A.S. Hedayat, N.J.A. Sloane and J. Stufken, 1999, Orthogonal Arrays: Theory and Applications, Springer Series in Statistics.

F. Gamboa, A. Janon, T. Klein and A. Lagnoux, 2014, Sensitivity indices for multivariate outputs, Electronic Journal of Statistics, 8:575-603.

J.Y. Tissot and C. Prieur, 2015, Estimating Sobol's indices combining Monte Carlo integration and Latin hypercube sampling, J. Statist. Comput. Simulation, 85:1358-1381.

See Also

sobolmara, sobolroauc, sobolMultOut

Examples

library(boot)
library(numbers)

####################
# Test case: the non-monotonic Sobol g-function

# The method of sobol requires 2 samples
# (there are 8 factors, all following the uniform distribution on [0,1])

# first-order sensitivity indices
x <- sobolroalhs(model = sobol.fun, factors = 8, N = 1000, order = 1, nboot=100)
print(x)
plot(x)

# closed second-order sensitivity indices
x <- sobolroalhs(model = sobol.fun, factors = 8, N = 1000, order = 2, nboot=100)
print(x)
plot(x)

####################
# Test case: dealing with non-uniform distributions

x <- sobolroalhs(model = NULL, factors = 3, N = 1000, order =1, nboot=0)

# X1 follows a log-normal distribution:
x$X[,1] <- qlnorm(x$X[,1])

# X2 follows a standard normal distribution:
x$X[,2] <- qnorm(x$X[,2])

# X3 follows a gamma distribution:
x$X[,3] <- qgamma(x$X[,3],shape=0.5)

# toy example
toy <- function(x){rowSums(x)}
y <- toy(x$X)
tell(x, y)
print(x)
plot(x)

####################
# Test case : multidimensional outputs

toy <- function(x){cbind(x[,1]+x[,2]+x[,1]*x[,2],2*x[,1]+3*x[,1]*x[,2]+x[,2])}
x <- sobolroalhs(model = toy, factors = 3, N = 1000, p=2, order =1, nboot=100)
print(x)
plot(x)