Search code examples
r

sensitivity - Why the response has more rows compared to the design


I want to use the sensitivity::tell(x, y) but the required y is longer than x

I made a Sobol design X and then (outside R) I executed the experiments to get the y. Now I have the design and the data. Looking at this response seems that I should use the sensitivity::tell but the length of y doesn't make sense.

The example code goes like that

n <- 1000
X1 <- data.frame(matrix(runif(8 * n), nrow = n))
X2 <- data.frame(matrix(runif(8 * n), nrow = n))

#Sensitivity analysis
x <- sobol(model = NULL, X1, X2, nboot = 100)
y <- ishigami.fun(x$X)
tell(x,y)
print(x)
plot(x)

It works, but the ishigami.fun(x$X) gives a vector much longer compared to the input design. If the design is of length 2,000 (X1 & X2) how the response could be 9,000 observations? If my design has 2000 experiments-rows the response should be equally long.

How can I get the sensitivity of my model using these functions?

enter image description here


Solution

  • The sobol function uses something called the pick-freeze estimator. It requires two designs (X1 and X2) because it is going to replace columns of one design with columns of the other in various ways, and it stacks these permuted designs (see "Calculation of indices" section of Wikipedia article). The final size of the stacked designs is the size of the original design (X1 or X2) times (number of effects + 1). Since you have 8 columns of X1 and the sobol function defaults to only first order effects, this results in a design with 1000 x (8+1) = 9000 rows.

    The trick to making this workflow work is to not run your function at your original design, but to run it after you obtain your pick-freeze estimator design x$X, and then use the tell function.

    If you have already done your function evaluations and you did not use a design that will work with the sobol function, you may want to consider using those runs to build an emulator (or surrogate model), and then get the Sobol decomposition of the emulator. For example, if you use the BASS package to build an emulator, you can use BASS::sobol to get the Sobol decomposition. For the Ishigami functions you could do this via

    n <- 1000
    X <- data.frame(matrix(runif(8 * n), nrow = n))
    y <- sensitivity::ishigami.fun(X)
    
    library(BASS)
    mod <- bass(X, y) # fit surrogate
    plot(mod) # check that surrogate fits well
    sob <- BASS::sobol(mod) # get Sobol decomposition of surrogate
    plot(sob)