qosa | R Documentation |
qosa
implements the estimation of first-order quantile-oriented sensitivity indices
as defined in Fort et al. (2016) with a kernel-based estimator of conditonal probability density functions
closely related to the one proposed by Maume-Deschamps and Niang (2018).
qosa
also supports a kernel-based estimation of Sobol first-order indices (i.e. Nadaraya-Watson).
qosa(model = NULL, X1, X2 = NULL, type = "quantile", alpha = 0.1, split.sample = 2/3,
nsample = 1e4, nboot = 0, conf = 0.95, ...)
## S3 method for class 'qosa'
tell(x, y = NULL, ...)
## S3 method for class 'qosa'
print(x, ...)
## S3 method for class 'qosa'
plot(x, ylim = c(0, 1), ...)
## S3 method for class 'qosa'
ggplot(data, mapping = aes(), ylim = c(0, 1), ..., environment
= parent.frame())
model |
a function, or a model with a |
X1 |
a random sample of the inputs used for the estimation of conditional probability density functions.
If |
X2 |
a random sample of the inputs used to evaluate the conditional probability density functions.
If NULL, it is constructed with the last |
type |
a string specifying which first-order sensitivity indices must be estimated: quantile-oriented indices ( |
alpha |
if |
split.sample |
if |
nsample |
the number of samples from the conditional probability density functions used
to estimate the conditional quantiles (if |
nboot |
the number of bootstrap replicates. |
conf |
the confidence level for confidence intervals. |
x |
a list of class |
data |
a list of class |
y |
a vector of model responses. |
ylim |
y-coordinate plotting limits. |
mapping |
Default list of aesthetic mappings to use for plot. If not specified, must be supplied in each layer added to the plot. |
environment |
[Deprecated] Used prior to tidy evaluation. |
... |
any other arguments for |
Quantile-oriented sensitivty indices were defined as a special case of sensitivity indices based on contrast functions in Fort et al. (2016).
The estimator used by qosa
follows closely the one proposed by Maume-Deschamps & Niang (2018).
The only difference is that Maume-Deschamps and Niang (2018) use the following kernel-based estimate of the conditional cumulative distribution function:
\hat{F}(y\Vert X=x) = \frac{ \sum_{i=1}^n K_{h_x}(x - X_i) \bold{1}\{Y_i< y\}}{\sum_{i=1}^n K_{h_x}(x - X_i)}
whereas we use
\hat{F}(y\vert X=x) = \frac{ \sum_{i=1}^n K_{h_x}(x - X_i) \int_{-\infty}^y K_{h_y}(t - Y_i)dt} {\sum_{i=1}^n K_{h_x}(x - X_i)},
meaning that \bold{1}\{Y_i< y\}
is replaced by \int_{-\infty}^y K_{h_y}(t - Y_i)dt = \Phi (\frac{y-Y_i}{h_y})
where \Phi
is the cumulative distribution function of the standard normal distribution (since kernel K
is Gaussian).
The two definitions thus coincide when h_y \rightarrow 0
. Our formula arises from a kernel density estimator of the joint pdf with a diagonal bandwidth.
In a future version, it will be genralized to a general bandwidth matrix for improved performance.
qosa
returns a list of class "qosa"
, containing all
the input arguments detailed before, plus the following components:
call |
the matched call. |
X |
a |
X1 |
a |
X |
a |
y |
a vector of model responses. |
S |
the estimations of the Sobol' sensitivity indices. |
Sebastien Da Veiga
Fort, J. C., Klein, T., and Rachdi, N. (2016). New sensitivity analysis subordinated to a contrast. Communications in Statistics-Theory and Methods, 45(15), 4349-4364.
Maume-Deschamps, V., and Niang, I. (2018). Estimation of quantile oriented sensitivity indices. Statistics & Probability Letters, 134, 122-127.
library(ks)
library(ggplot2)
library(boot)
# Test case : difference of two exponential distributions (Fort et al. (2016))
# We use two samples with different sizes
n1 <- 5000
X1 <- data.frame(matrix(rexp(2 * n1,1), nrow = n1))
n2 <- 1000
X2 <- data.frame(matrix(rexp(2 * n2,1), nrow = n2))
Y1 <- X1[,1] - X1[,2]
Y2 <- X2[,1] - X2[,2]
x <- qosa(model = NULL, X1, X2, type = "quantile", alpha = 0.1)
tell(x,c(Y1,Y2))
print(x)
ggplot(x)
# Test case : difference of two exponential distributions (Fort et al. (2016))
# We use only one sample
n <- 1000 # put n=10000 for more consistency
X <- data.frame(matrix(rexp(2 * n,1), nrow = n))
Y <- X[,1] - X[,2]
x <- qosa(model = NULL, X1 = X, type = "quantile", alpha = 0.7)
tell(x,Y)
print(x)
ggplot(x)
# Test case : the Ishigami function
# We estimate first-order Sobol' indices (by specifying 'mean')
# Next lines are put in comment because too long fro CRAN tests
#n <- 5000
#nboot <- 50
#X <- data.frame(matrix(-pi+2*pi*runif(3 * n), nrow = n))
#x <- qosa(model = ishigami.fun, X1 = X, type = "mean", nboot = nboot)
#print(x)
#ggplot(x)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.