sensitivity (version 1.10.1)

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

Description

sobolroalhs implements the estimation of the Sobol' sensitivity indices introduced by Tissot & Prieur (2012) using two Orthogonal Array-based Latin Hypercubes. This function allows the estimation of all first-order indices or all closed second-order indices (containing the sum of the second-order effect between two inputs and the individual effects of each input) at a total cost of $2 \times N$. For closed second-order indices $N=q^{2}$ where $q \geq d-1$ is a prime number denoting the number of levels of the orthogonal array, and where $d$ is the number of factors.

Usage

sobolroalhs(model = NULL, factors, runs, order, conf=0.95, tail=TRUE, na.rm=FALSE, ...)
## S3 method for class 'sobolroalhs':
tell(x, y = NULL, \dots)
## S3 method for class 'sobolroalhs':
print(x, \dots)
## S3 method for class 'sobolroalhs':
plot(x, ylim = c(0,1), type="standard", ...)

Arguments

model

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

factors

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

runs

an integer specifying the number $N$ of model runs.

order

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

conf

the confidence level for confidence intervals.

tail

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

na.rm

a boolean specifying if the response of the model contains NA values.

x

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

y

a vector of model responses.

ylim

coordinate plotting limits for the indices values.

type

a character specifying the type of estimator to plot (standard for the basic estimator or monod for the Janon-Monod estimator.)

...

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:
• call
  the matched call.
• X
  a matrix containing the design of experiments.
• OA
  the orthogonal array constructed (NULL if order=1).
• levels
  the number of levels of the orthogonal array constructed (NULL if order=1).
• y
  a vector of model responses.
• V
  a data.frame containing the estimations of the variance (Vs for the standard variance and Veff for the Janon-Monod variance).
• S
  a data.frame containing the estimations of the Sobol' sensitivity indices (S for the standard estimator and Seff for the Janon-Monod estimator).

Details

The method used by sobolroalhs only considers models whose inputs follow uniform distributions on [0,1]. The transformations of each input (between its initial distribution and a U[0,1] distribution) have therefore to be realized before the call to sobolroalhs(). Bootstrap confidence intervals are not available with this method ; the given confidence intervals come from asymptotical formula.

References

Tissot, J. Y. and Prieur, C. (2012), Estimating Sobol's indices combining Monte Carlo integration and Latin hypercube sampling http://hal.archives-ouvertes.fr/hal-00743964. A.S. Hedayat, N.J.A. Sloane, John Stufken (1999), Orthogonal Arrays: Theory and Applications.

Examples

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, runs = 1000, order = 1)
print(x)
plot(x)

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

# Test case : the Ishigami function

# New function because sobolroalhs() works with U[0,1] inputs
ishigami1.fun=function(x) ishigami.fun(x*2*pi-pi)

# first-order sensitivity indices
x <- sobolroalhs(model = ishigami1.fun, factors = 3, runs = 100000, order = 1)
print(x)
plot(x)

# closed second-order sensitivity indices
x <- sobolroalhs(model = ishigami1.fun, factors = 3, runs = 100000, order = 2)
print(x)
plot(x)

# dealing with NA values
x <- sobolroalhs(model = NULL, factors = 3, runs = 100000, order =1,na.rm=TRUE)
y <- ishigami1.fun(x$X)
# we randomly insert NA values in y
pos <- sample(length(y),100)
y[pos] <- NA
tell(x,y)
print(x)
plot(x)