mixle.doe.sensitivity module

Global sensitivity analysis: variance-based Sobol indices and Morris screening.

Which inputs actually drive a model’s output? Sobol indices decompose the output variance into the contribution of each input (first order) and of each input including all its interactions (total order). Morris screening is a cheaper one-at-a-time elementary-effects method for an initial factor ranking. These tell you which survey parameters / forcings to refine and which to fix – the front of the UQ loop.

sobol_indices(func, bounds, n=4096, *, seed=0, names=None, n_bootstrap=200, confidence=0.95)[source]

First- and total-order Sobol sensitivity indices (Saltelli sampling, Jansen estimators).

Parameters:

  • func (Callable[[ndarray], ndarray]) – a vectorized model f(X) -> y mapping an (m, d) array of inputs to an (m,) array of scalar outputs.

  • bounds (Sequence[tuple[float, float]]) – [(lo, hi), ...] per input – the input is taken uniform on the box.

  • n (int) – base sample size; the total number of model evaluations is n * (d + 2) (bootstrap resampling below is free of additional func calls, so this count is exact).

  • seed (int) – RNG seed for the (Sobol) base samples.

  • names (Sequence[str] | None) – optional input names for the returned dict.

  • n_bootstrap (int) – number of bootstrap resamples used for the standard error / confidence interval of each index (see “Returns” below).

  • confidence (float) – confidence level for the bootstrap interval, in (0, 1) – e.g. 0.95 for a two-sided 95% interval.

Returns:

  • 'S1' / 'ST': the raw, unclipped point estimate. A finite-sample Sobol estimator can legitimately land slightly negative, or (for S1) slightly above 1, purely from Monte Carlo noise – most often when the true index is near 0 or 1. That is not a bug in the estimate; it is the honest signal that sampling noise is comparable to the index’s true magnitude, and this raw value is always what is returned here (never silently clamped).

  • 'S1_clipped' / 'ST_clipped': the same estimate clipped into its theoretically valid range (S1 into [0, 1], ST below at 0) – a separate, explicitly-named convenience for e.g. plotting or reporting a single non-negative “share of variance”; never the only value returned, and never used in place of the raw estimate above.

  • 'S1_standard_error' / 'ST_standard_error': a bootstrap standard error, from n_bootstrap resamples (with replacement) of the already-evaluated A/B/AB_i rows, each re-run through the same S1/ST formula – the standard uncertainty-quantification approach for Sobol’ estimators (Archer, Saltelli & Sobol’ 1997), since no exact closed-form sampling distribution exists for them in general.

  • 'S1_ci_low' / 'S1_ci_high' and 'ST_ci_low' / 'ST_ci_high': the corresponding percentile bootstrap confidence interval at the confidence level.

Plus 'names', 'var' (the total output variance), and 'sampling_method' (normally 'scrambled_sobol'; 'pseudorandom_fallback' only when SciPy’s QMC capability is absent). S1[i] is the fraction of output variance from input i alone; ST[i] includes all interactions involving i (so ST[i] - S1[i] measures i’s interaction strength, and ST[i] ~ 0 – relative to its standard error, not just numerically close – means input i can be fixed).

Return type:

A dict with, for both first- and total-order

morris_screening(func, bounds, *, trajectories=20, levels=4, seed=0, names=None)[source]

Morris elementary-effects screening – a low-cost first factor ranking.

Walks trajectories randomized-direction one-factor-at-a-time paths on a levels-point grid per dimension (Morris 1991; Campolongo, Cariboni & Saltelli 2007); the mean absolute elementary effect mu_star[i] ranks influence and the spread sigma[i] flags nonlinearity/interactions. Cost: trajectories * (d + 1) evaluations – far fewer than Sobol, for an initial screen.

Each trajectory walks the grid {0, 1/(p-1), 2/(p-1), ..., 1} (p = levels) and takes d one-at-a-time steps, one per input in a random order, each changing exactly one coordinate by the fixed step delta = p / (2*(p-1)) = (p/2) / (p-1) – i.e. exactly p/2 grid points, chosen (the standard choice) so the step always lands exactly on another grid point. Both the visit order and each dimension’s step direction (up by delta or down by delta) are independently randomized per trajectory. The starting coordinate for a dimension stepping up is drawn only from the grid’s lower half (indices 0..p/2-1, which by construction have room for a +delta move without leaving [0, 1]); for a dimension stepping down, only from the upper half (indices p/2..p-1, symmetric room for -delta). Every step therefore lands on the grid and inside bounds by construction – never clipped or silently shrunk to a smaller step – and, because the two halves are equal-sized and a dimension’s pre-step value is confined to its assigned half while its post-step value always lands in the other, every grid level is visited with equal (1/p) probability, marginally, by any one coordinate at any position along a trajectory: the design is unbiased across the input box, not concentrated away from either boundary.

levels must be even, so the grid-index step levels // 2 is exact (an odd level count has no step that both equals the standard delta formula and lands on a grid point).

Parameters:
Return type:

dict[str, Any]

fast_indices(func, bounds, n=600, *, harmonics=6, seed=0, names=None, n_bootstrap=200, confidence=0.95)[source]

First-order sensitivity indices via Random Balance Designs FAST (RBD-FAST).

A Fourier alternative to sobol_indices() for the first-order indices: every input is driven along the same triangle-wave search curve but under an independent random permutation, the model is evaluated once over the n points, and for each input the output – reordered along that input’s curve – has its variance concentrated at the base frequency’s first harmonics harmonics. The ratio of that power to the total is the first-order index (with the Tarantola bias correction). Cost is a single batch of n evaluations, independent of dimension.

S1 is the raw Tarantola-corrected finite-sample estimate. S1_clipped is a separate presentation convenience. Circular moving-block bootstrap fields report the raw estimator’s standard error and percentile interval without additional model evaluations.

Parameters:
Return type:

dict[str, Any]

dgsm(func, bounds, n=1024, *, seed=0, rel_step=1.0e-4, names=None)[source]

Derivative-based global sensitivity measures (DGSM): mean squared partial derivatives.

nu[i] = E[(df/dx_i)^2] over the input box, estimated by central finite differences at n low-discrepancy points. Unlike the first-order Sobol index, a DGSM is nonzero whenever an input matters anywhere – including purely through interactions – so it is a low-cost, robust screen that upper-bounds the total Sobol index (Sobol & Kucherenko, via the Poincare inequality: ST[i] <= (L_i / pi)^2 * nu[i] / Var(y) for a uniform input of width L_i). The reported importance is L_i^2 * nu[i] normalized to sum to one – a dimensionless influence ranking.

Returns {'nu': (d,), 'importance': (d,), 'names': [...]}.

Parameters:
Return type:

dict[str, Any]