Optimal Transport Global Sensitivity Analysis

• The topic presents a framework that leverages Wasserstein distances to measure shifts in full output distributions, replacing traditional scalar summaries.
• It unifies methodologies across PDEs, stochastic computer codes, and causal interventions by embedding distributional outputs in Wasserstein space.
• Computational strategies include LP solvers and rank-based estimators for practical implementation, addressing challenges like cubic complexity and estimation bias.

Searching arXiv for the specified papers and closely related work on optimal transport-based sensitivity analysis. Optimal Transport-based Global Sensitivity Analysis denotes a family of sensitivity-analysis frameworks in which discrepancies between outputs generated under different inputs, parameter values, or intervention rules are quantified by optimal-transport functionals, most commonly Wasserstein distances. Across the literature, the approach appears in several distinct but related settings: sensitivity of PDE solutions such as the LWR traffic-flow model on networks (Briani et al., 2016); global sensitivity analysis for computer codes whose outputs are probability measures or cumulative distribution functions, and for stochastic computer codes (Fort et al., 2020); sensitivity analysis for Fokker–Planck equations through quantitative Wasserstein bounds (Morange, 3 Feb 2026); sensitivity analysis for stochastic interventions in causal inference via Monge–Kantorovich formulations (Levis et al., 2024); and model-agnostic post-processing tools for OT-based sensitivity indices implemented in software (Chiani et al., 24 Jul 2025). In all of these formulations, the central idea is to replace scalar dispersion summaries by a transport-based comparison of full output laws, thereby measuring “mass shift” rather than only pointwise or moment-wise error (Briani et al., 2016).

1. Conceptual scope and problem classes

Optimal Transport-based Global Sensitivity Analysis is not a single index but a collection of methodologies adapted to the geometry of distribution-valued, stochastic, or law-dependent outputs. In the formulation for black-box codes, one considers a code

$Z = f(X_1,\ldots,X_p)$

with independent inputs $X_1,\ldots,X_p$ and output $Z$ taking values in a metric space $(\mathcal X,d)$; a principal special case is $\mathcal X=\mathcal W_q(\mathbb R)$, the space of real probability measures with finite $q$th moment, endowed with the $q$-Wasserstein distance (Fort et al., 2020). In this setting, sensitivity indices are defined directly on the space of output laws rather than on scalar summaries.

A second problem class concerns stochastic computer codes. There the code has the form

$f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$

where $d$ is an unobserved “random seed.” For each fixed $x$, the law of $X_1,\ldots,X_p$0 is a probability measure $X_1,\ldots,X_p$1, which induces the “ideal deterministic” code

$X_1,\ldots,X_p$2

(Fort et al., 2020). This construction embeds stochastic simulation outputs into Wasserstein space and permits a standard global-sensitivity analysis at the level of output distributions.

A third class is PDE- and transport-model sensitivity. For traffic networks, the LWR density on each directed edge satisfies

$X_1,\ldots,X_p$3

with junction dynamics encoded through distribution matrices $X_1,\ldots,X_p$4 and sensitivity measured by Wasserstein distance between density fields generated by different inputs (Briani et al., 2016). For diffusion-type models, the Fokker–Planck equation

$X_1,\ldots,X_p$5

is analyzed with respect to parameter perturbations using explicit upper bounds on $X_1,\ldots,X_p$6 (Morange, 3 Feb 2026).

A fourth class arises in causal sensitivity analysis under stochastic or generalized interventions. With observed data $X_1,\ldots,X_p$7, generalized policies induce couplings between observational and target treatment distributions, and worst-case sensitivity bounds can be written as Monge–Kantorovich problems over couplings $X_1,\ldots,X_p$8 (Levis et al., 2024). This suggests that OT-based sensitivity analysis is equally a framework for uncertainty propagation and for partial-identification geometry.

2. Mathematical formulations in Wasserstein space

The unifying object is an optimal-transport cost between probability measures. In the Kantorovich formulation, for measures $X_1,\ldots,X_p$9 on $Z$0 and nonnegative cost $Z$1,

$Z$2

When $Z$3, $Z$4 is the $Z$5-Wasserstein cost (Chiani et al., 24 Jul 2025). In one dimension,

$Z$6

for c.d.f.s $Z$7 (Fort et al., 2020).

For OT-based sensitivity indices in model-agnostic GSA, if $Z$8 is the model output and $Z$9 an input, a “separation” functional $(\mathcal X,d)$0 satisfies $(\mathcal X,d)$1, and the local separation is

$(\mathcal X,d)$2

The associated global index is

$(\mathcal X,d)$3

Using OT as the separation, one sets

$(\mathcal X,d)$4

with numerator

$(\mathcal X,d)$5

Normalization uses

$(\mathcal X,d)$6

where $(\mathcal X,d)$7 is an independent copy, and yields the first-order OT sensitivity index

$(\mathcal X,d)$8

(Chiani et al., 24 Jul 2025).

In the Wasserstein-space Sobol-style construction, if $(\mathcal X,d)$9 is a random element in $\mathcal X=\mathcal W_q(\mathbb R)$0 with Fréchet mean

$\mathcal X=\mathcal W_q(\mathbb R)$1

then the Fréchet variance is

$\mathcal X=\mathcal W_q(\mathbb R)$2

For a subset $\mathcal X=\mathcal W_q(\mathbb R)$3, the conditional Fréchet mean $\mathcal X=\mathcal W_q(\mathbb R)$4 defines the sensitivity index

$\mathcal X=\mathcal W_q(\mathbb R)$5

In one dimension for $\mathcal X=\mathcal W_q(\mathbb R)$6,

$\mathcal X=\mathcal W_q(\mathbb R)$7

(Fort et al., 2020).

For networked conservation laws, if $\mathcal X=\mathcal W_q(\mathbb R)$8 are nonnegative densities on a network $\mathcal X=\mathcal W_q(\mathbb R)$9 with equal total mass $q$0, then

$q$1

where $q$2 is the network geodesic distance (Briani et al., 2016). For treatment-policy sensitivity, the worst-case deviation from the identified component $q$3 is

$q$4

and minimizing worst-case width over $q$5 yields the Monge–Kantorovich problem

$q$6

with dual

$q$7

(Levis et al., 2024).

3. Principal index constructions and structural properties

Several OT-based indices have been proposed, reflecting different output objects and inferential goals. For distribution-valued outputs, the Wasserstein–Fréchet approach uses the Fréchet mean and Fréchet variance in $q$8 and yields Sobol-style indices $q$9 with $q$0, $q$1, and invariance under translation, isometries or nondegenerate scaling of the output (Fort et al., 2020). A second construction in the same work is based on Hoeffding decomposition of indicator functions of Wasserstein balls. For fixed $q$2 and target distributions $q$3,

$q$4

Integrating the associated partial variances over reference pairs leads to

$q$5

where

$q$6

(Fort et al., 2020).

The model-agnostic OT index implemented in software is based on conditional-vs-marginal distributional separation. Its normalization is designed so that $q$7 and $q$8 almost surely, with monotonicity and normalization properties (Chiani et al., 24 Jul 2025). This framework is especially adapted to multivariate outputs $q$9 and does not require independence of the inputs $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$0 (Chiani et al., 24 Jul 2025).

When the cost is squared Euclidean, the index admits a decomposition using the Gelbrich formula:

$f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$1

leading to

$f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$2

The Wasserstein-Bures semi-metric index is $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$3 (Chiani et al., 24 Jul 2025). This decomposition connects OT-based sensitivity to mean effects, covariance effects, and higher-order moment effects.

In PDE settings, the role of the sensitivity index is often played by a normalized transport discrepancy rather than a variance ratio. For the LWR model, given a baseline input vector $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$4 and perturbed input $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$5, one computes

$f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$6

where $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$7 is the discrete transport optimum, and records either $f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$8 or the time-average

$f_s: E\times D\to\mathbb R,\quad (x,d)\mapsto y,$9

(Briani et al., 2016). For Fokker–Planck equations, Morange’s framework defines a transport-based sensitivity index

$d$0

relative to a nominal $d$1, and proves $d$2 (Morange, 3 Feb 2026).

4. Computational procedures and algorithms

A core practical feature of OT-based sensitivity analysis is the reduction of sensitivity computation to optimal-transport solvers. For the LWR network model, each edge is discretized into uniform cells $d$3 of length $d$4, producing an undirected graph $d$5 whose nodes correspond to cell centers. Masses are accumulated as

$d$6

with equal total mass $d$7, and the cost matrix is built from shortest-path distances $d$8, for example by Dijkstra. The resulting Hitchcock linear program is

$d$9

subject to

$x$0

Then $x$1, with error bound $x$2 (Briani et al., 2016).

For generic input-output samples, the implemented OT-GSA workflow partitions the domain of an input $x$3 into $x$4 bins and, for each bin, computes an empirical OT cost between the marginal output sample and the conditional output sample restricted to that bin. The numerator is averaged over bins, and the result is normalized by an empirical upper bound

$x$5

to obtain $x$6 (Chiani et al., 24 Jul 2025). The special cases of one-dimensional outputs and Gaussian/Bures geometry admit closed forms handled by ot_indices_1d and ot_indices_wb (Chiani et al., 24 Jul 2025).

Two solvers are emphasized in software. Network Simplex solves the classical OT linear program with worst-case complexity roughly $x$7, while Sinkhorn–Knopp solves the entropic-regularized problem

$x$8

with per-iteration cost $x$9 and linear convergence in the sense that $X_1,\ldots,X_p$00 exponentially fast; the “sinkhorn-stable” variant avoids underflow/overflow for small $X_1,\ldots,X_p$01, and $X_1,\ldots,X_p$02 as $X_1,\ldots,X_p$03 (Chiani et al., 24 Jul 2025).

For stochastic computer codes in Wasserstein spaces, practical estimation proceeds by generating $X_1,\ldots,X_p$04 outer samples of the inputs and, for each outer sample, $X_1,\ldots,X_p$05 inner runs corresponding to different random seeds to form empirical measures

$X_1,\ldots,X_p$06

These empirical measures are then plugged into pick-freeze estimators, U-statistic estimators, or rank-based formulas for $X_1,\ldots,X_p$07 or $X_1,\ldots,X_p$08 (Fort et al., 2020).

For Fokker–Planck sensitivity, implementation is analytic rather than solver-centered. One first estimates or bounds $X_1,\ldots,X_p$09, $X_1,\ldots,X_p$10, and the uniform ellipticity $X_1,\ldots,X_p$11, then assembles

$X_1,\ldots,X_p$12

and uses the resulting bound to control output-law perturbations over a parameter-uncertainty radius (Morange, 3 Feb 2026).

5. Quantitative bounds, analytical results, and canonical couplings

The Fokker–Planck setting provides an explicit global-in-time Lipschitz-in-parameter transport estimate. Under the uniform-ellipticity and global-Lipschitz assumptions, for all $X_1,\ldots,X_p$13,

$X_1,\ldots,X_p$14

If the initial law is unchanged, this becomes

$X_1,\ldots,X_p$15

(Morange, 3 Feb 2026). Two proof strategies yield the same constants: synchronous coupling of the underlying SDEs and differentiation of the Kantorovich dual (Morange, 3 Feb 2026).

In the overdamped Langevin case,

$X_1,\ldots,X_p$16

with $X_1,\ldots,X_p$17 and $X_1,\ldots,X_p$18, one has the sharper estimate

$X_1,\ldots,X_p$19

where $X_1,\ldots,X_p$20 (Morange, 3 Feb 2026). The separation of contraction and parameter-shift effects is a distinctive feature of this result.

In causal sensitivity analysis, the sharp worst-case bounds depend on the cost structure. Under Model 2, where

$X_1,\ldots,X_p$21

the optimal cost is

$X_1,\ldots,X_p$22

so that

$X_1,\ldots,X_p$23

As $X_1,\ldots,X_p$24, the bounds collapse to $X_1,\ldots,X_p$25 (Levis et al., 2024). Under Model 3, where $X_1,\ldots,X_p$26, the minimal cost is the $X_1,\ldots,X_p$27-Wasserstein distance,

$X_1,\ldots,X_p$28

and

$X_1,\ldots,X_p$29

(Levis et al., 2024).

The optimal couplings in these causal models are analytically characterized. For Model 2, the optimal coupling is the maximal coupling $X_1,\ldots,X_p$30 with

$X_1,\ldots,X_p$31

implemented by a two-stage policy that keeps $X_1,\ldots,X_p$32 with probability $X_1,\ldots,X_p$33 and otherwise redraws $X_1,\ldots,X_p$34 (Levis et al., 2024). For Model 3 with continuous $X_1,\ldots,X_p$35 and $X_1,\ldots,X_p$36, the optimal transport for convex $X_1,\ldots,X_p$37 is the increasing rearrangement

$X_1,\ldots,X_p$38

(Levis et al., 2024). These results show that OT-based sensitivity analysis can identify not only bound widths but also extremal generalized policies.

6. Applications, numerical findings, and limitations

The LWR traffic-flow study is a canonical application of OT-based sensitivity analysis to large networks. The baseline input vector $X_1,\ldots,X_p$39 collects the initial condition $X_1,\ldots,X_p$40, fundamental-diagram parameters $X_1,\ldots,X_p$41, the distribution matrices $X_1,\ldots,X_p$42, and the network topology $X_1,\ldots,X_p$43; perturbations are introduced factor by factor, and the transport discrepancy $X_1,\ldots,X_p$44 is evaluated over time (Briani et al., 2016). The reported numerical findings are specific: initial-condition errors produced by small-scale shifts yield a $X_1,\ldots,X_p$45 that decays rapidly in time, so $X_1,\ldots,X_p$46 is small for large $X_1,\ldots,X_p$47; sensitivity to $X_1,\ldots,X_p$48 and $X_1,\ldots,X_p$49 grows roughly linearly with $X_1,\ldots,X_p$50 and $X_1,\ldots,X_p$51, with a proportionality constant increasing with network size; a small bias in a junction distribution coefficient causes significant $X_1,\ldots,X_p$52 growth in time and saturation; perturbing all junctions amplifies sensitivity roughly proportionally to network diameter; and topology changes such as closing one central road produce sensitivity of order one, with long-time $X_1,\ldots,X_p$53 scaling linearly with network size (Briani et al., 2016). The paper further states that the LP-based computation is robust but scales as $X_1,\ldots,X_p$54 in the worst case, so moderate grid resolution and network size are recommended in practice (Briani et al., 2016).

For distribution-valued and stochastic computer codes, numerical studies illustrate both first-level and second-level GSA. The toy cdf model permits closed-form expressions for $X_1,\ldots,X_p$55 and $X_1,\ldots,X_p$56, while empirical comparisons show that rank-based estimators often have substantially smaller MSE than pick-freeze estimators at equal code-call budgets (Fort et al., 2020). In the stochastic-code version of the same toy model, empirical output measures are built from i.i.d. draws and the rank-based method again outperforms pick-freeze for first-order indices (Fort et al., 2020). Second-level GSA is formulated by treating the distributions of the inputs themselves as uncertain and computing OT-based indices with respect to input-law parameters (Fort et al., 2020).

The software package gsaot operationalizes the model-agnostic OT-GSA framework as a post-processing step requiring only an input matrix x and an output matrix y; it supports ot_indices, ot_indices_wb, ot_indices_1d, ot_indices_smap, irrelevance_threshold, and associated plotting methods (Chiani et al., 24 Jul 2025). Practical examples include a linear Gaussian test, a spruce budworm ODE system, and a climate module with a custom $X_1,\ldots,X_p$57 cost matrix (Chiani et al., 24 Jul 2025). The package emphasizes multivariate outputs, correlated inputs, and visualization of local separation curves (Chiani et al., 24 Jul 2025).

Several limitations recur across the literature. In the Fokker–Planck theory, the current results cover $X_1,\ldots,X_p$58; the $X_1,\ldots,X_p$59 case might be handled via reflection coupling, while the $X_1,\ldots,X_p$60 case remains open (Morange, 3 Feb 2026). The drift and diffusion must be globally Lipschitz, and $X_1,\ldots,X_p$61 uniformly elliptic; singular or degenerate diffusions require new ideas (Morange, 3 Feb 2026). In generic OT solvers, worst-case complexity is cubic, entropic regularization introduces bias for finite $X_1,\ldots,X_p$62, and binning smooths local structure (Chiani et al., 24 Jul 2025). Rank-based methods for Wasserstein-space GSA are not yet available for higher-order indices or vector-valued inputs, and U-statistic kernels can become high-dimensional (Fort et al., 2020).

A common misconception is that transport-based methods merely repackage variance-based sensitivity. The cited works indicate otherwise. OT-based indices are formulated on full conditional and marginal laws, can be applied to multivariate outputs and correlated inputs, and in the causal setting yield sharp nonparametric bound widths that vanish as the target policy approaches the observational regime (Levis et al., 2024, Chiani et al., 24 Jul 2025). A plausible implication is that optimal transport serves not only as a discrepancy metric but as a structural language for global sensitivity whenever the output object is itself distributional, stochastic, or intervention-dependent.