gsaot gsaot website

R-CMD-check test-coverage CRAN status

The package gsaot provides a set of tools to compute and plot Optimal Transport (OT) based sensitivity indices. The core functions of the package are:

The package also provides functions to plot the resulting indices and the separation measures.

Installation

install.packages("gsaot")

You can install the development version of gsaot from GitHub with:

# install.packages("remotes")
remotes::install_github("pietrocipolla/gsaot")

:exclamation: :exclamation: Installation note

The sinkhorn and sinkhorn_log solvers in gsaot greatly benefit from optimization in compilation. To add this option (before package installation), edit your .R/Makevars file with the desired flags. Even though different compilers have different options, a common flag to enable a safe level of optimization is

CXXFLAGS+=-O2

More detailed information on how to customize the R packages compilation can be found in the R guide.

Example

We can use a gaussian toy model with three outputs as an example:

library(gsaot)

N <- 1000

mx <- c(1, 1, 1)
Sigmax <- matrix(data = c(1, 0.5, 0.5, 0.5, 1, 0.5, 0.5, 0.5, 1), nrow = 3)

x1 <- rnorm(N)
x2 <- rnorm(N)
x3 <- rnorm(N)

x <- cbind(x1, x2, x3)
x <- mx + x %*% chol(Sigmax)

A <- matrix(data = c(4, -2, 1, 2, 5, -1), nrow = 2, byrow = TRUE)
y <- t(A %*% t(x))

x <- data.frame(x)

After having defined the number of partitions, we compute the sensitivity indices using different solvers. First, Sinkhorn solver and default parameters:

M <- 25

sensitivity_indices <- ot_indices(x, y, M)
sensitivity_indices
#> Method: sinkhorn 
#> 
#> Indices:
#>        X1        X2        X3 
#> 0.6040598 0.6304113 0.2817856

Second, Network Simplex solver:

sensitivity_indices <- ot_indices(x, y, M, solver = "transport")
sensitivity_indices
#> Method: transport 
#> 
#> Indices:
#>        X1        X2        X3 
#> 0.5135492 0.5274847 0.1734944

Third, Wasserstein-Bures solver, with bootstrap:

sensitivity_indices <- ot_indices_wb(x, y, M, boot = TRUE, R = 100)
sensitivity_indices
#> Method: wass-bures 
#> 
#> Indices:
#>        X1        X2        X3 
#> 0.4833691 0.4940299 0.1097725 
#> 
#> Advective component:
#>        X1        X2        X3 
#> 0.2943062 0.3174669 0.1018953 
#> 
#> Diffusive component:
#>          X1          X2          X3 
#> 0.189062870 0.176562994 0.007877195 
#> 
#> Type of confidence interval: norm 
#> Number of replicates: 100 
#> Confidence level: 0.95 
#> Bootstrap statistics:
#>   input  component   original        bias      low.ci    high.ci
#> 1    X1 wass-bures 0.49284730 0.009478185 0.465069799 0.50166843
#> 2    X2 wass-bures 0.50378617 0.009756254 0.479159970 0.50889985
#> 3    X3 wass-bures 0.12863517 0.018862723 0.087497514 0.13204739
#> 4    X1  advective 0.29925643 0.004950185 0.282197530 0.30641496
#> 5    X2  advective 0.32163503 0.004168111 0.308057464 0.32687637
#> 6    X3  advective 0.11239018 0.010494927 0.083363925 0.12042659
#> 7    X1  diffusive 0.19359087 0.004528000 0.181197928 0.19692781
#> 8    X2  diffusive 0.18215114 0.005588143 0.169509166 0.18361682
#> 9    X3  diffusive 0.01624499 0.008367797 0.002543453 0.01321094

Fourth, we can use the package to compute the sensitivity map on the output:

sensitivity_indices <- ot_indices_smap(x, y, M)
sensitivity_indices
#>             X1         X2        X3
#> [1,] 0.5865925 0.04945865 0.2004293
#> [2,] 0.3221281 0.71433091 0.1206923

mirror server hosted at Truenetwork, Russian Federation.