Last data update: 2014.03.03

R: Sobol' Indices Estimation Using Replicated OA-based LHS
sobolroalhsR Documentation

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*N. For closed second-order indices N=q^2 where q >= 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, ...)
## S3 method for class 'sobolroalhs'
print(x, ...)
## 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.

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.

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).

Warning messages

"The number of model evaluations (runs) you entered is not the square of a prime number. It has been replaced by : "

when order=2, the number of levels of the orthogonal array must be a prime number. If the number of runs specified is not a square of a prime number then this warning message indicates that the number of runs was replaced depending on the value of tail. If tail=TRUE (resp. tail=FALSE) the new number of runs is equals to the square of the prime number preceding (resp. following) the square root of runs.

"The number of model evaluations (runs) you entered is not satisfying the constraint n >= (d-1)^2. It has been replaced by : "

when order=2, the number of runs must satisfied the constraint N ≥ (d-1)^2 where d is the number of factors. This warning message indicates that the number of runs was replaced by the square of the prime number following (or equals to) d-1.

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.

See Also

sobolmara

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)

Results