### Install GlobalSensitivity.jl

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/index.md

Install the GlobalSensitivity.jl package using the Julia package manager.

```julia
using Pkg
Pkg.add("GlobalSensitivity")
```

--------------------------------

### Sobol Method Examples

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Demonstrates using the Sobol method with parameter ranges, explicit design matrices, and batched evaluation. Access first-order (S1), second-order (S2), and total-order (ST) indices.

```julia
using GlobalSensitivity, QuasiMonteCarlo

function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

lb = -ones(3) * π
ub =  ones(3) * π

# Method 1: using parameter ranges directly
res1 = gsa(ishi, Sobol(order = [0, 1, 2]), [[lb[i], ub[i]] for i in 1:3], samples = 10000)
println("S1 = ", res1.S1)   # first-order indices
println("S2 = ", res1.S2)   # second-order interaction indices
println("ST = ", res1.ST)   # total-order indices

# Method 2: using explicit design matrices (more control over sampling)
samples = 100000
A, B = QuasiMonteCarlo.generate_design_matrices(samples, lb, ub, SobolSample())
res2 = gsa(ishi, Sobol(nboot = 10, conf_level = 0.95), A, B)
println("S1 confidence intervals = ", res2.S1_Conf_Int)

# Method 3: batched evaluation
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end
res3 = gsa(ishi_batch, Sobol(), A, B, batch = true)
```

--------------------------------

### View Julia Version Information

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/index.md

Displays detailed version information about the Julia installation, including the operating system and architecture. Useful for reproducibility.

```julia
using InteractiveUtils # hide
versioninfo() # hide
```

--------------------------------

### eFAST Method Example

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Estimates first-order and total-order sensitivity indices using the eFAST method. Supports uniform bounds, arbitrary input distributions, and batched computations.

```julia
using GlobalSensitivity, Distributions

function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

lb = -ones(3) * π; ub = ones(3) * π

# uniform bounds
res1 = gsa(ishi, eFAST(), [[lb[i], ub[i]] for i in 1:3], samples = 15000)
println("S1 = ", res1.S1)
println("ST = ", res1.ST)

# with arbitrary input distributions
input_dists = [Normal(0, 1), Uniform(-π, π), Uniform(-π, π)]
res2 = gsa(ishi, eFAST(), input_dists, samples = 15000)

# batched
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end
res3 = gsa(ishi_batch, eFAST(), [[lb[i], ub[i]] for i in 1:3], samples = 15000, batch = true)
```

--------------------------------

### DGSM Method Example

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Computes derivative-based global sensitivity measures using DGSM. Supports arbitrary input distributions and can compute cross-derivatives when `crossed=true`.

```julia
using GlobalSensitivity, Distributions

f1(x) = x[1] + 2 * x[2] + 6.0 * x[3]

# each parameter can have its own distribution
dist1 = [Uniform(4, 10), Normal(4, 23), Beta(2, 3)]

res = gsa(f1, DGSM(), dist1, samples = 100_000)
println("Mean gradient:       ", res.a)
println("Mean |gradient|:     ", res.absa)
println("Mean gradient²:      ", res.asq)
println("Sigma:               ", res.sigma)
println("Tao:                 ", res.tao)

# with cross-derivatives (uses Hessian via ForwardDiff)
res_crossed = gsa(f1, DGSM(crossed = true), dist1, samples = 10_000)
println("Crossed indices:\n", res_crossed.crossed)
```

--------------------------------

### DeltaMoment Method Example

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Computes moment-independent density-based sensitivity indices using DeltaMoment. Supports scalar output and can be used with pre-sampled data.

```julia
using GlobalSensitivity, QuasiMonteCarlo

function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

lb = -ones(3) * π; ub = ones(3) * π
p_range = fill([lb[1], ub[1]], 3)

res = gsa(ishi, DeltaMoment(), p_range, samples = 1000)
println("Delta indices:          ", res.deltas)
println("Adjusted delta indices: ", res.adjusted_deltas)
println("CI lower:               ", res.adjusted_deltas_low)
println("CI upper:               ", res.adjusted_deltas_hi)

# from pre-sampled data
X = QuasiMonteCarlo.sample(1000, lb, ub, QuasiMonteCarlo.SobolSample())
Y = [ishi(X[:, i]) for i in axes(X, 2)]
res2 = gsa(X, Y, DeltaMoment())
```

--------------------------------

### Setup Uncorrelated Parameter Distributions

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/shapley.md

Defines the marginal distributions and a Gaussian Copula with an identity covariance matrix to represent uncorrelated input parameters for sensitivity analysis.

```julia
d = length(θ)
mu = zeros(Float32, d)
#covariance matrix for the copula
Covmat = Matrix(1.0f0 * I, d, d)
#the marginal distributions for each parameter
marginals = [Normal(mu[i]) for i in 1:d]

copula = GaussianCopula(Covmat)
input_distribution = SklarDist(copula, marginals)
```

--------------------------------

### Setup Correlated Parameter Distributions

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/shapley.md

Defines a Gaussian Copula with a specified correlation matrix (0.09 off-diagonal) to represent parameters with interdependencies for sensitivity analysis.

```julia
Corrmat = fill(0.09f0, d, d)
for i in 1:d
    Corrmat[i, i] = 1.0f0
end

#since the marginals are standard normal the covariance matrix and correlation matrix are the same
copula = GaussianCopula(Corrmat)
input_distribution = SklarDist(copula, marginals)
```

--------------------------------

### Lotka-Volterra Model Setup and Regression GSA

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Sets up the Lotka-Volterra ODE problem and performs Global Sensitivity Analysis using the Regression GSA method. Visualizes the results using heatmaps for partial and standard regression coefficients.

```julia
using GlobalSensitivity, QuasiMonteCarlo, OrdinaryDiffEq, Statistics, CairoMakie

function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2] #prey
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2] #predator
end

u0 = [1.0; 1.0]
tspan = (0.0, 10.0)
p = [1.5, 1.0, 3.0, 1.0]
prob = ODEProblem(f, u0, tspan, p)
t = collect(range(0, stop = 10, length = 200))

f1 = function (p)
    prob1 = remake(prob; p = p)
    sol = solve(prob1, Tsit5(); saveat = t)
    return [mean(sol[1, :]), maximum(sol[2, :])]
end

bounds = [[1.0, 5.0], [1.0, 5.0], [1.0, 5.0], [1.0, 5.0]]

reg_sens = gsa(f1, RegressionGSA(true), bounds, samples = 200)
fig = Figure(resolution = (600, 400))
ax,
hm = CairoMakie.heatmap(fig[1, 1],
    reg_sens.partial_correlation,
    axis = (xticksvisible = false, yticksvisible = false, yticklabelsvisible = false,
        xticklabelsvisible = false, title = "Partial correlation"))
Colorbar(fig[1, 2], hm)
ax,
hm = CairoMakie.heatmap(fig[2, 1],
    reg_sens.standard_regression,
    axis = (xticksvisible = false, yticksvisible = false, yticklabelsvisible = false,
        xticklabelsvisible = false, title = "Standard regression"))
Colorbar(fig[2, 2], hm)
fig

```

--------------------------------

### Morris Method Example with ODEProblem

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Applies the Morris method to a Lotka-Volterra predator-prey model defined as an ODEProblem. It calculates elementary effects on the mean and maximum population sizes.

```julia
using GlobalSensitivity, OrdinaryDiffEq, Statistics, StableRNGs

function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2]   # prey
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2]  # predator
end
u0 = [1.0; 1.0]; tspan = (0.0, 10.0); p = [1.5, 1.0, 3.0, 1.0]
prob = ODEProblem(f, u0, tspan, p)
t = collect(range(0, stop = 10, length = 200))

f1 = function (p)
    sol = solve(remake(prob; p = p), Tsit5(); saveat = t)
    [mean(sol[1, :]), maximum(sol[2, :])]
end

_rng = StableRNG(1234)
m = gsa(f1, Morris(total_num_trajectory = 1000, num_trajectory = 150),
        [[1.0, 5.0], [1.0, 5.0], [1.0, 5.0], [1.0, 5.0]], rng = _rng)
```

--------------------------------

### RegressionGSA Method Example

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Computes linear and rank-based sensitivity measures using the RegressionGSA method. Requires defining the model function, parameter bounds, and optionally providing pre-sampled matrices.

```julia
using GlobalSensitivity, QuasiMonteCarlo, OrdinaryDiffEq, Statistics

function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2]
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2]
end
u0 = [1.0; 1.0]; prob = ODEProblem(f, u0, (0.0, 10.0), [1.5, 1.0, 3.0, 1.0])
t = collect(range(0, 10, length = 200))

f1 = p -> (sol = solve(remake(prob; p = p), Tsit5(); saveat = t);
           [mean(sol[1, :]), maximum(sol[2, :])])

bounds = [[1.0, 5.0], [1.0, 5.0], [1.0, 5.0], [1.0, 5.0]]

# with rank-based coefficients
reg = gsa(f1, RegressionGSA(true), bounds, samples = 200)
println("Pearson:              ", reg.pearson)
println("SRC:                  ", reg.standard_regression)
println("PCC:                  ", reg.partial_correlation)
println("Spearman (rank):      ", reg.pearson_rank)

# from pre-sampled matrices
lb = [1.0, 1.0, 1.0, 1.0]; ub = [5.0, 5.0, 5.0, 5.0]
X = QuasiMonteCarlo.sample(500, lb, ub, QuasiMonteCarlo.SobolSample())
Y = reduce(hcat, [f1(X[:, i]) for i in axes(X, 2)])  # 2×500 output matrix
reg_mat = gsa(X, Y, RegressionGSA(true))
```

--------------------------------

### View Package Status

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/index.md

Displays the status of installed Julia packages. This is often used for reproducibility checks.

```julia
using Pkg # hide
Pkg.status() # hide
```

--------------------------------

### Import necessary packages

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Import the required libraries for Global Sensitivity Analysis, statistics, ODE solving, quasi-Monte Carlo sampling, and plotting.

```julia
using GlobalSensitivity, Statistics, OrdinaryDiffEq, QuasiMonteCarlo, Plots
```

--------------------------------

### Initialize RSA Method

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Initializes the RSA method. Use `n_dummy_parameters` to set the number of dummy parameters and `acceptance_threshold` to define the criterion for classifying outputs.

```julia
RSA(; n_dummy_parameters::Int = 10, acceptance_threshold = mean)
```

--------------------------------

### GSA with Design Matrices and Batch Execution

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Demonstrates using pre-generated design matrices for Sobol GSA and performing analysis in a batched mode for performance. Includes timing comparisons between standard and batched execution.

```julia
using QuasiMonteCarlo
samples = 500
lb = [1.0, 1.0, 1.0, 1.0]
ub = [5.0, 5.0, 5.0, 5.0]
sampler = SobolSample()
A, B = QuasiMonteCarlo.generate_design_matrices(samples, lb, ub, sampler)
sobol_sens_desmat = gsa(f1, Sobol(), A, B)

f_batch = function (p)
    prob_func(prob, i, repeat) = remake(prob; p = p[:, i])
    output_func(sol, i) = ([mean(sol[1, :]), maximum(sol[2, :])], false)
    ensemble_prob = EnsembleProblem(prob; prob_func, output_func)

    sol = solve(
        ensemble_prob, Tsit5(), EnsembleThreads(); saveat = t, trajectories = size(p, 2))

    out = reshape(sol, :, size(p, 2))

    return out
end

sobol_sens_batch = gsa(f_batch, Sobol(), A, B, batch = true)

@time gsa(f1, Sobol(), A, B)
@time gsa(f_batch, Sobol(), A, B, batch = true)

```

--------------------------------

### Define the Lotka-Volterra ODE model

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Define the Lotka-Volterra equations, initial conditions, time span, and parameters to set up the ODE problem.

```julia
function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2] #prey
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2] #predator
end
u0 = [1.0; 1.0]
tspan = (0.0, 10.0)
p = [1.5, 1.0, 3.0, 1.0]
prob = ODEProblem(f, u0, tspan, p)
t = collect(range(0, stop = 10, length = 200))
```

--------------------------------

### Generate Sobol design matrices with QuasiMonteCarlo

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Manually generate the design matrices 'A' and 'B' for the Sobol method using `QuasiMonteCarlo.jl`. This allows for greater control over the sampling process and the use of specific quasi-random sequences.

```julia
samples = 500
lb = [1.0, 1.0, 1.0, 1.0]
ub = [5.0, 5.0, 5.0, 5.0]
sampler = SobolSample()
A, B = QuasiMonteCarlo.generate_design_matrices(samples, lb, ub, sampler)
```

--------------------------------

### Sobol Method: Ishigami Function (Serial)

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/README.md

Performs Sobol sensitivity analysis on the Ishigami function using serial execution. Requires defining the function, bounds, and sampler.

```julia
function ishi(X)
    A = 7
    B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

n = 600000
lb = -ones(3) * π
ub = ones(3) * π
sampler = SobolSample()
A, B = QuasiMonteCarlo.generate_design_matrices(n, lb, ub, sampler)

res1 = gsa(ishi, Sobol(order = [0, 1, 2]), A, B)
```

--------------------------------

### View Manifest and Project Status

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/index.md

Shows the status of all dependencies, including those in the manifest file, for a comprehensive view of the project's environment. This is useful for detailed reproducibility.

```julia
using Pkg # hide
Pkg.status(; mode = PKGMODE_MANIFEST) # hide
```

--------------------------------

### Plotting First Order Sobol Indices

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Visualizes first-order Sobol indices for 'Prey' and 'Predator' using scatter plots. Requires CairoMakie and pre-computed Sobol indices.

```julia
ax,
hm = CairoMakie.scatter(
    fig[1, 2], sobol_sens.S1[2][1, 2:end], label = "Prey", markersize = 4)
CairoMakie.scatter!(
    fig[1, 2], sobol_sens.S1[2][2, 2:end], label = "Predator", markersize = 4)
```

```julia
ax,
hm = CairoMakie.scatter(
    fig[2, 1], sobol_sens.S1[3][1, 2:end], label = "Prey", markersize = 4)
CairoMakie.scatter!(
    fig[2, 1], sobol_sens.S1[3][2, 2:end], label = "Predator", markersize = 4)
```

```julia
ax,
hm = CairoMakie.scatter(
    fig[2, 2], sobol_sens.S1[4][1, 2:end], label = "Prey", markersize = 4)
CairoMakie.scatter!(
    fig[2, 2], sobol_sens.S1[4][2, 2:end], label = "Predator", markersize = 4)
```

--------------------------------

### Sobol GSA with Bootstrapping for Confidence Intervals

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Performs Sobol Global Sensitivity Analysis with bootstrapping to estimate confidence intervals for the sensitivity indices. Visualizes the first-order and total-order indices for prey and predator.

```julia
f1 = function (p)
    prob1 = remake(prob; p = p)
    sol = solve(prob1, Tsit5(); saveat = t)
end
sobol_sens = gsa(f1, Sobol(nboot = 20), bounds, samples = 500)
fig = Figure(resolution = (600, 400))
ax,
hm = CairoMakie.scatter(
    fig[1, 1], sobol_sens.S1[1][1, 2:end], label = "Prey", markersize = 4)
CairoMakie.scatter!(fig[1, 1], sobol_sens.S1[1][2, 2:end], label = "Predator", markersize = 4)


```

--------------------------------

### Sobol Method: Ishigami Function (Batch)

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/README.md

Performs Sobol sensitivity analysis on the Ishigami function using the batching interface for potentially faster execution. The function must be vectorized.

```julia
function ishi_batch(X)
    A = 7
    B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end

res2 = gsa(ishi_batch, Sobol(), A, B, batch = true)
```

--------------------------------

### General GSA Interface

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

The unified entry point for all GSA methods. `f` maps parameters to output, `param_range` defines parameter bounds, and `samples` controls evaluation count. Set `batch=true` for matrix-based input/output.

```julia
gsa(f, method::GSAMethod, param_range; samples, batch = false)
```

```julia
using GlobalSensitivity, QuasiMonteCarlo

# scalar output function
function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

# batched version (operates on columns)
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end

lb = -ones(3) * π
ub =  ones(3) * π
p_range = [[lb[i], ub[i]] for i in 1:3]

# Serial call
res_serial = gsa(ishi, Sobol(), p_range, samples = 10000)

# Batched (parallel) call
res_batch = gsa(ishi_batch, Sobol(), p_range, samples = 10000, batch = true)

println(res_serial.S1)   # First-order Sobol indices
println(res_serial.ST)   # Total-order Sobol indices
```

--------------------------------

### Initialize Mutual Information Method

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Initializes the Mutual Information method. `n_bootstraps` controls the number of bootstrap samples for the null distribution, and `conf_level` sets the confidence level for the bounds.

```julia
MutualInformation(; n_bootstraps = 1000, conf_level = 0.95)
```

--------------------------------

### General Interface - gsa

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

The unified entry point for all GSA methods. It takes a function `f`, a GSA method, parameter ranges, and optional arguments like `samples` and `batch` for parallel evaluation.

```APIDOC
## General Interface - gsa

### Description
The unified entry point for all methods.

### Signature
```julia
gsa(f, method::GSAMethod, param_range; samples, batch = false)
```

### Parameters
- `f`: A function mapping a parameter vector to a scalar or vector output.
- `method`: An instance of a `GSAMethod` (e.g., `Sobol()`, `Morris()`).
- `param_range`: A vector of `[lb, ub]` pairs defining the range for each parameter.
- `samples`: Controls the number of evaluations.
- `batch`: When `true`, `f` receives a parameter matrix and must return a matrix.

### Example
```julia
using GlobalSensitivity, QuasiMonteCarlo

# scalar output function
function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

# batched version (operates on columns)
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end

lb = -ones(3) * π
ub =  ones(3) * π
p_range = [[lb[i], ub[i]] for i in 1:3]

# Serial call
res_serial = gsa(ishi, Sobol(), p_range, samples = 10000)

# Batched (parallel) call
res_batch = gsa(ishi_batch, Sobol(), p_range, samples = 10000, batch = true)

println(res_serial.S1)   # First-order Sobol indices
println(res_serial.ST)   # Total-order Sobol indices
```
```

--------------------------------

### Sobol Method Configuration

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Configure the Sobol method to compute specific orders of sensitivity indices and bootstrap for confidence intervals. `Ei_estimator` selects the Monte Carlo estimator for total indices.

```julia
Sobol(; order = [0, 1], nboot = 1, conf_level = 0.95)
```

--------------------------------

### Perform GSA with Parallel ODE Solves using GlobalSensitivity.jl

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/joss/paper.md

This snippet demonstrates how to perform Global Sensitivity Analysis (GSA) on an ODE problem using GlobalSensitivity.jl. It leverages the Ensemble Interface with EnsembleThreads for multithreaded parallelization of ODE solves and calculates the mean of prey and maximum of predator populations as outputs. The Sobol sampler is used for generating design matrices for the GSA.

```julia
using GlobalSensitivity, QuasiMonteCarlo, OrdinaryDiffEq

function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2] #prey
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2] #predator
end

u0 = [1.0; 1.0]
tspan = (0.0, 10.0)
p = [1.5, 1.0, 3.0, 1.0]
prob = ODEProblem(f, u0, tspan, p)
t = collect(range(0, stop = 10, length = 200))

f1 = function (p)
    prob_func(prob, i, repeat) = remake(prob; p = p[:, i])
    ensemble_prob = EnsembleProblem(prob, prob_func = prob_func)
    sol = solve(
        ensemble_prob, Tsit5(),
        EnsembleThreads();
        saveat = t, trajectories = size(p, 2))
    # Now sol[i] is the solution for the ith set of parameters
    out = zeros(2, size(p, 2))
    for i in 1:size(p, 2)
        out[1, i] = mean(sol[i][1, :])
        out[2, i] = maximum(sol[i][2, :])
    end
    out
end

samples = 10000
lb = [1.0, 1.0, 1.0, 1.0]
ub = [5.0, 5.0, 5.0, 5.0]
sampler = SobolSample()
A, B = QuasiMonteCarlo.generate_design_matrices(samples, lb, ub, sampler)

sobol_result = gsa(f1, Sobol(), A, B, batch = true)
```

--------------------------------

### Morris Method Configuration

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Configuration options for the Morris method, a screening method for parameter importance.

```APIDOC
## Morris Method

### Description
Screening method (OAT — One At a Time) that computes elementary effects to rank parameters by importance. Returns `means`, `means_star` (absolute means), and `variances` of the elementary effects. High `means_star` indicates importance; high variance indicates nonlinearity or interaction.

### Constructor
```julia
Morris(; p_steps = Int[], relative_scale = false, num_trajectory = 10,
         total_num_trajectory = 5 * num_trajectory, len_design_mat = 10)
```

### Parameters
- `p_steps`: Number of steps for each parameter in the design matrix.
- `relative_scale`: Whether to use relative scaling for parameters.
- `num_trajectory`: Number of trajectories to generate.
- `total_num_trajectory`: Total number of trajectories (defaults to 5 * `num_trajectory`).
- `len_design_mat`: Length of the design matrix.

### Example
```julia
using GlobalSensitivity, OrdinaryDiffEq, Statistics, StableRNGs

function f(du, u, p, t)
    du[1] = p[1] * u[1] - p[2] * u[1] * u[2]   # prey
    du[2] = -p[3] * u[2] + p[4] * u[1] * u[2]  # predator
end
u0 = [1.0; 1.0]; tspan = (0.0, 10.0); p = [1.5, 1.0, 3.0, 1.0]
prob = ODEProblem(f, u0, tspan, p)
t = collect(range(0, stop = 10, length = 200))

f1 = function (p)
    sol = solve(remake(prob; p = p), Tsit5(); saveat = t)
    [mean(sol[1, :]), maximum(sol[2, :])]
end

_rng = StableRNG(1234)
m = gsa(f1, Morris(total_num_trajectory = 1000, num_trajectory = 150),
        [[1.0, 5.0], [1.0, 5.0], [1.0, 5.0], [1.0, 5.0]], rng = _rng)
```
```

--------------------------------

### Sobol and eFAST GSA Visualization for Prey Population

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Visualizes the first-order and total-order sensitivity indices for the prey population using Sobol and eFAST methods. Bar plots are used to display the results.

```julia
fig = Figure(resolution = (600, 400))
barplot(fig[1, 1],
    [1, 2, 3, 4],
    sobol_sens.S1[1, :],
    color = :green,
    axis = (xticksvisible = false, xticklabelsvisible = false,
        title = "Prey (Sobol)", ylabel = "First order"))
barplot(fig[2, 1], [1, 2, 3, 4], sobol_sens.ST[1, :], color = :green,
    axis = (xticksvisible = false, xticklabelsvisible = false, ylabel = "Total order"))
barplot(fig[1, 2], [1, 2, 3, 4], efast_sens.S1[1, :], color = :red,
    axis = (xticksvisible = false, xticklabelsvisible = false, title = "Prey (eFAST)"))
barplot(fig[2, 2], [1, 2, 3, 4], efast_sens.ST[1, :], color = :red,
    axis = (xticksvisible = false, xticklabelsvisible = false))
fig

```

--------------------------------

### Perform RSA with Custom Threshold

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Performs Global Sensitivity Analysis using the RSA method with a custom acceptance threshold set to the 75th percentile of the output.

```julia
using GlobalSensitivity

function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end

lb = -ones(3) * π; ub = ones(3) * π

res1 = gsa(ishi, RSA(), [[lb[i], ub[i]] for i in 1:3], samples = 15000)
println("Sensitivity indices: ", res1.S)
println("Dummy param stats (mean, std): ", res1.Sd)

# custom acceptance threshold (e.g., 75th percentile)
res2 = gsa(ishi, RSA(acceptance_threshold = Y -> quantile(Y, 0.75)),
           [[lb[i], ub[i]] for i in 1:3], samples = 15000)

# batched
res3 = gsa(ishi_batch, RSA(), [[lb[i], ub[i]] for i in 1:3], samples = 15000, batch = true)
```

--------------------------------

### Perform Sobol Global Sensitivity Analysis

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Conduct a Sobol sensitivity analysis on the function `f1` with specified parameter ranges and number of samples. This method provides first-order and total-order effects.

```julia
m = gsa(f1, Sobol(), [[1.0, 5.0], [1.0, 5.0], [1.0, 5.0], [1.0, 5.0]], samples = 1000)
```

--------------------------------

### Morris Method Configuration

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Configure the Morris method for screening parameters. `num_trajectory` and `total_num_trajectory` control the number of trajectories, influencing the granularity of elementary effect calculations.

```julia
Morris(; p_steps = Int[], relative_scale = false, num_trajectory = 10,
         total_num_trajectory = 5 * num_trajectory, len_design_mat = 10)
```

--------------------------------

### Generate Dataset for ODE

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/shapley.md

Generates simulated data for an Ordinary Differential Equation (ODE) model. This involves defining the ODE function, time span, and solving the problem to obtain the data.

```julia
using GlobalSensitivity, OrdinaryDiffEq, Flux, SciMLSensitivity, LinearAlgebra
using Optimization, OptimizationOptimisers, Distributions, Copulas, CairoMakie

u0 = [2.0f0; 0.0f0]
datasize = 30
tspan = (0.0f0, 1.5f0)

function trueODEfunc(du, u, p, t)
    true_A = [-0.1f0 2.0f0; -2.0f0 -0.1f0]
    du .= ((u .^ 3)'true_A)'
end
t = range(tspan[1], tspan[2], length = datasize)
prob = ODEProblem(trueODEfunc, u0, tspan)
ode_data = Array(solve(prob, Tsit5(), saveat = t))
```

--------------------------------

### Sobol and eFAST GSA Visualization for Predator Population

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/juliacon21.md

Visualizes the first-order and total-order sensitivity indices for the predator population using Sobol and eFAST methods. Bar plots are used to display the results.

```julia
fig = Figure(resolution = (600, 400))
barplot(fig[1, 1],
    [1, 2, 3, 4],
    sobol_sens.S1[2, :],
    color = :green,
    axis = (xticksvisible = false, xticklabelsvisible = false,
        title = "Predator (Sobol)", ylabel = "First order"))
barplot(fig[2, 1], [1, 2, 3, 4], sobol_sens.ST[2, :], color = :green,
    axis = (xticksvisible = false, xticklabelsvisible = false, ylabel = "Total order"))
barplot(fig[1, 2], [1, 2, 3, 4], efast_sens.S1[2, :], color = :red,
    axis = (xticksvisible = false, xticklabelsvisible = false, title = "Predator (eFAST)"))
barplot(fig[2, 2], [1, 2, 3, 4], efast_sens.ST[2, :], color = :red,
    axis = (xticksvisible = false, xticklabelsvisible = false))
fig

```

--------------------------------

### Define Optimization Problem

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/shapley.md

Sets up an optimization problem using Optimization.jl to solve a cost function derived from an ODE. This is a prerequisite for sensitivity analysis.

```julia
adtype = Optimization.AutoZygote()

optf = Optimization.OptimizationFunction((p, _) -> loss_n_ode(p), adtype)
optprob = Optimization.OptimizationProblem(optf, θ)

result_neuralode = Optimization.solve(optprob,
    OptimizationOptimisers.Adam(0.05),
    callback = callback,
    maxiters = 300)
```

--------------------------------

### Sobol Method Configuration

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

Configuration options for the Sobol method, a variance-based GSA technique.

```APIDOC
## Sobol Method

### Description
Variance-based method that decomposes output variance into contributions from individual parameters and their interactions, returning first-order (`S1`), second-order (`S2`), and total-order (`ST`) indices.

### Constructor
```julia
Sobol(; order = [0, 1], nboot = 1, conf_level = 0.95)
```

### Parameters
- `order`: Indices to compute — `[0,1]` for first+total, add `2` for second-order interactions.
- `nboot`: Bootstrap runs for confidence intervals (requires `nboot > 1`).
- `conf_level`: Confidence level for confidence intervals.

### Example
```julia
using GlobalSensitivity, QuasiMonteCarlo

function ishi(X)
    A = 7; B = 0.1
    sin(X[1]) + A * sin(X[2])^2 + B * X[3]^4 * sin(X[1])
end

lb = -ones(3) * π
ub =  ones(3) * π

# Method 1: using parameter ranges directly
res1 = gsa(ishi, Sobol(order = [0, 1, 2]), [[lb[i], ub[i]] for i in 1:3], samples = 10000)
println("S1 = ", res1.S1)   # first-order indices
println("S2 = ", res1.S2)   # second-order interaction indices
println("ST = ", res1.ST)   # total-order indices

# Method 2: using explicit design matrices (more control over sampling)
samples = 100000
A, B = QuasiMonteCarlo.generate_design_matrices(samples, lb, ub, SobolSample())
res2 = gsa(ishi, Sobol(nboot = 10, conf_level = 0.95), A, B)
println("S1 confidence intervals = ", res2.S1_Conf_Int)

# Method 3: batched evaluation
function ishi_batch(X)
    A = 7; B = 0.1
    @. sin(X[1, :]) + A * sin(X[2, :])^2 + B * X[3, :]^4 * sin(X[1, :])
end
res3 = gsa(ishi_batch, Sobol(), A, B, batch = true)
```
```

--------------------------------

### FractionalFactorial Method for Main Effects

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

This method uses a Hadamard-matrix design to efficiently compute main effects and their variances. The design size is determined by the number of parameters, and no explicit `samples` argument is needed. It can be used with explicit parameter ranges.

```julia
FractionalFactorial()
```

```julia
using GlobalSensitivity

# interaction between params 7 and 12 only
f = X -> X[1] + 2 * X[2] + 3 * X[3] + 4 * X[7] * X[12]

main_effects, variances = gsa(f, FractionalFactorial(), num_params = 12)
println("Main effects: ", main_effects)
println("Variances:    ", variances)

# with explicit parameter ranges
p_range = [[-1.0, 1.0] for _ in 1:12]
main_effects2, _ = gsa(f, FractionalFactorial(), num_params = 12, p_range = p_range)
```

--------------------------------

### Create a function for sensitivity analysis

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Define a function that takes parameters, solves the ODE, and returns the average prey population and maximum predator population. This function is used as the model for GSA.

```julia
f1 = function (p)
    prob1 = remake(prob; p = p)
    sol = solve(prob1, Tsit5(); saveat = t)
    [mean(sol[1, :]), maximum(sol[2, :])]
end
```

--------------------------------

### Plot Sobol indices

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Visualize the total order (ST) and first-order (S1) Sobol indices for both the average prey and maximum predator populations. This helps in understanding the contribution of each parameter to the model output variance.

```julia
p1 = bar(["a", "b", "c", "d"], sobol_result.ST[1, :],
    title = "Total Order Indices prey", legend = false)
p2 = bar(["a", "b", "c", "d"], sobol_result.S1[1, :],
    title = "First Order Indices prey", legend = false)
p1_ = bar(["a", "b", "c", "d"], sobol_result.ST[2, :],
    title = "Total Order Indices predator", legend = false)
p2_ = bar(["a", "b", "c", "d"], sobol_result.S1[2, :],
    title = "First Order Indices predator", legend = false)
plot(p1, p2, p1_, p2_)
```

--------------------------------

### RBDFAST Method for First-Order Indices

Source: https://context7.com/sciml/globalsensitivity.jl/llms.txt

RBDFAST is an efficient first-order sensitivity index estimator using Random Balance Designs and random permutations. It supports both serial and batched function evaluations. Ensure to provide a `StableRNG` for reproducibility.

```julia
RBDFAST(; num_harmonics = 6)
```

```julia
using GlobalSensitivity, StableRNGs

function linear(X)
    A = 7; B = 0.1
    A * X[1] + B * X[2]
end
function linear_batch(X)
    A = 7; B = 0.1
    @. A * X[1, :] + B * X[2, :]
end

rng = StableRNG(123)
res1 = gsa(linear, GlobalSensitivity.RBDFAST(), num_params = 2, samples = 15000, rng = rng)
println("First-order indices (serial):  ", res1)

rng2 = StableRNG(123)
res2 = gsa(linear_batch, GlobalSensitivity.RBDFAST(), num_params = 2,
           samples = 15000, batch = true, rng = rng2)
println("First-order indices (batched): ", res2)
```

--------------------------------

### Perform Parallelized GSA with Sobol Method

Source: https://github.com/sciml/globalsensitivity.jl/blob/master/docs/src/tutorials/parallelized_gsa.md

Executes the global sensitivity analysis using the Sobol method in parallel by passing `batch = true` to the `gsa` function. This enables multithreaded computation of the sensitivity indices.

```julia
sobol_result = gsa(f1, Sobol(), A, B, batch = true)
```