Counterfactual Identifiability via Dynamic Optimal Transport
Authors: Fabio De Sousa Ribeiro
Ainkaran Santhirasekaram
Ben Glocker
Venue: NeurIPS 2025 poster ()
Keywords: generative models
causality
counterfactual inference
optimal transport
continuous normalizing flows
flow matching
Primary area: deep_learning
TL;DR: We address the open question of counterfactual identification for high-dimensional multivariate outcomes from observational data.

Abstract: We address the open question of counterfactual identification for high-dimensional multivariate outcomes from observational data. Pearl (2000) argues that counterfactuals must be identifiable (i.e., recoverable from the observed data distribution) to justify causal claims. A recent line of work on counterfactual inference shows promising results but lacks identification, undermining the causal validity of its estimates. To address this, we establish a foundation for multivariate counterfactual identification using continuous-time flows, including non-Markovian settings under standard criteria. We characterise the conditions under which flow matching yields a unique, monotone and rank-preserving counterfactual transport map with tools from dynamic optimal transport, ensuring consistent inference. Building on this, we validate the theory in controlled scenarios with counterfactual ground-truth and demonstrate improvements in axiomatic counterfactual soundness on real images.

Counterfactual Identifiability via Dynamic Optimal Transport
Fabio De Sousa Ribeiro
Ainkaran Santhirasekaram
Ben Glocker
Imperial College London, UK
Abstract
We address the open question of counterfactual identification for high-dimensional multivariate outcomes from observational data. Pearl (2000) argues that counter- factuals must be identifiable (i.e., recoverable from the observed data distribution) to justify causal claims. A recent line of work on counterfactual inference shows promising results but lacks identification, undermining the causal validity of its estimates. To address this, we establish foundation for multivariate counterfactual identification using continuous-time flows, including non- Mon-Markovian settings under standard criteria. We characterise the conditions under which flow matching yields a unique, monotone, and rank touk-preeervng counterfactual transport map with tools from dynamic optimal transport, ensuring consistent inference. Building on this, we validate the theory in controlled scenarios with counterfactual ground truth and demonstrate improvements in axiomatic counterfactual soundness real images.
Introduction
has been argued that shortcomings today's deep models reveal an overreliance on statistical associations and lack of causal understanding (Schölkopf et al., 2021; Bareinboim al., 2022). This view has sparked interest in causality within machine learning (Peters et al., 2017: Castro et al., 2020; Schölkopf, 2022). Causal reasoning is the process of drawing conclusions from a causal model, which represents our assumptions about the data-generating process. Thus, causality and generative modelling are inextricably linked. Pearl (2009) formalised causal inference using the Structural Causal Model (SCM) framework, in which variables are generated by functional assignments and dependencies represented by graph. This framework (with do-calculus (Pearl, 1995)) expresses causal queries and, under stated assumptions, identifies causal effects that can be estimated from data.
Why is identification essential? Without identification it is impossible to make a precise causal claim: there can exist observationally equivalent models that yield different answers, and providing such guarantees is what causal methods are designed for (Bareinboim et al., 2022; Hyvärinen et al., 2024). Indeed, do-calculus (Pearl, 1995) is an identification tool that tells us whether a causal effect can be determined from observed data and set of assumptions encoded in causal graph. Further, Pearl (2000) insists that any counterfactual query be bound by an identifiability requirement, that is, that there exists unique mapping from the observed data distribution to the counterfactual of interest.
This work focuses on the identification of -ihh-imennsinall counterfactuals from observational data. Counterfactuals are hypothetical scenarios given observed evidence, for example, one may query "What would Y have been, had X been x. Counterfactuals have broad scientific utility, supporting the evaluation of interventions. (Kusner et al., 2017; Tsirtsis and Rodriguez, 2023), the characterisation of causal relationships (Karimi et al., 2020: Budhathoki et al., 2022), and the generation of targeted synthetic data for downstream tasks (Pitis et al., 2022; Roschewitz et al., 2024; Mehta et al., 2025).
Email: f.de-sousa-rrixiiro@@ addesoonssribbsiii@@mppaaal....
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

Deep generative models are increasingly being used to parameterise SCMs (Pawlowski et al., 2020; Sanchez and Tsaftaris, 2021; Ribeiro et al., 2023; Komanduri et al., 2024; Kumar al., 2025; Rasal et al., 2025; Xia et al., 2025). While the inferred counterfactuals are promising, the lack of identification precludes warranted causal interpretations of the estimates beyond statistical predictions. We aim to course-correct this practice. Counterfactual identification is arguably the most ambitious goal in causal analysis, as counterfactuals sit at the top of the causal hierarchy and subsume all other causal queries (Bareinboim et al., 2022). Thus, existing identification results are limited, with classical symbolic methods (Tian and Pearl, 2000; Pearl, 2001; Shpitser and Pearl, 2007) not having been developed for high-dimensiona variables (Xia et al., 2023; Nasr- Masr-Esfahany and Kiciman, 2023).
Nasr-Est Basr-Esfahany et al. (2023) established counterfactual identification of bijective mechanisms, propos- ing a flow-based estimation method similar to Pawlowski et al. (2020); Khemakhem et al. (2021). However, their identification results only extend to multi-dimensional variables under the Backdoor Criterion (BC) and sufficient variability (Hyvarinen et al., 2019), stating that: "it is not clear how to generalise the monotonicity condition to multi-ddimnnsioall variables" in Markovian causal struc- tures. This is well-known issue with close connections to ICA (Hyvärinen and Pajunen, 1999) and disentanglement (Locatello et al., 2019). In this work, we characterise a multivariate generalisation of monotonicity to provably solve this problem, without imposing arbitrary coordinate order.
We establish foundation for counterfactual identification of high-dimensional multivariate outcomes (e.g. images) from observational data, built on continuous time flows. We characterise the conditions under which flow matching yields a unique, rank- pre-erree counterfactual transport map using tools from dynamic optimal transport. Our results address causal validity issues in prior work and extend to Markovinn settings under standard criteria. We validate our theory in controlled scenarios with counterfactual ground-truth, and improve axiomatic counterfactual soundness on real images.
2 Related Work
Pearl's Causal Hierarchy (PCH) (Pearl and Mackenzie, 2018) defines three levels of causal abstraction: associational (C,); interventional (C2), and counterfactual (Cz). The PCH Theorem (Bareinboim et al., 2022) states that cross-layer inference of L: quantities using Cz data is generally not possible without further assumptions. Pearl (1995)'s do-calculus is sound and complete procedure for deter- mining the identifiability of an C2 causal quantity from Ln data, given a causal graph. Counterfactual identification (Lx) represents the most ambitious goal in causal amalysis as it sits at the top of the PCH and therefore subsumes all other causal queries (Pearl, 2009). Subsequently, known counterfactual identification results are comparatively limited in number and scope to date, and existing constraints needed for La identification can be broadly categorised as either graphical or functional in nature.
Classical identification techniques involve symbolic methods and assumptions about the causal graph to identify counterfactual queries from both L1 and L2 data (Tian and Pearl, 2000; Pearl, 2001; Shpitser and Pearl, 2007). Recently, graph-based criteria have been extended to nested counterfactuals and fairness-aware queries (Zhang and Bareinboim. 2018; Correa et al., 2021), whereas functional identification has moved towards neural and nonparametric schemes that approximate or bound counterfactuals (Hartford et al., 2017: Xia al., 2023; Melnychuk et al., 2023; Geffner et al., 2024; Pan and Bareinboim, 2024; Wu et al., 2025), even in non- non- -oniieentee settings (Gresele et al., 2022). Identification results predicated on monotonicity (Lu et al., 2020; Vlontzos et al., 2023; Javaloy et al., 2023) and bijectivity (Nasr-Esfahany et al., 2023) constraints are of primary interest to our work, since they generalise a broad range of model classes with known identifiability results, such as linear and nonlinear additive noise models (Shimizu et al., 2006; Hoyer et al., 2008; Peters et al., 2014), post-nonlinear models (Zhang and Hyvärinen, 2009), and location-scaa scale models (Immer et al., 2023).
Establishing counterfactual identification of high high-dimensional Markovian SCMs is a timely con- tribution, as many recent works rely on such assumptions (Pawlowski et al., 2020; Sanchez and Tsaftaris, 2021; Sánchez-Martin et al., 2022; Ribeiro et al., 2023; Chen et al., 2024; Rasal et al., 2025). However, these methods lack an identification strategy and cannot support causal claims, posing operational risk. Of particular relevance to our work is optimal transport (Brenier, 1991; Caffarelli, 1992; Benamou and Brenier, 2000; Villani et al., 2008; Santambrogio, 2015; Peyré et al., 2019), and modern continuous-time flow models (Chen et al., 2018; Lipman et al., 2023; Liu et al., 2023; Albergo and Vanden-Fin -ijnden,, 2023; Pooladian et al., 2023; Tong et al., 2024), which in combination form the basis of our theoretical results and practical counterfactual inference model prescription.
2

---

3 Preliminaries
3.1 Structural Causal Models
Definition 3.1 (Structural Causal Model (SCM) (Pearl, 2009)). An SCM e (U,X,F) is a
mathematical tool designed to express and infer causal quantities. It consists of: (i) a set of exogenous variables = {UI..... , Un} determined by factors outside of C. and distributed according to Pu: (ii) set of endogenous variables X {X1,...; Xn) determined by other variables in E; and (iii) set of functions F (fi..... fn} specifying the causal generative process mapping U to X:
X f(PAA, U.),
PA, X (X,),
for i = I,....
(I)
where PA, is the subset of endogenous variables that directly cause X,, called its parents. If the causal generative process is acyclic, an SCM can be represented by a Directed Acyclic Graph (DAG).
Counterfactual Inference. fact were different, all else being equal. SCMs can express and answer counterfactual queries of the form: "Given that we observed X x, what would X have been had X; been set to x" An intervention sets chosen variables to specified values, e.g. do(X; == x*) or dX == x(), where S (1...., n) indexes the intervened variables in the SCM, and the components of x's may differ'. Counterfactual inference proceeds in three steps: (i) Abduction: infer the posterior distribution over the exogenous variables PuIx=x- given observed evidence X = x; (ii) Action: intervene, e.g., apply do(X, == x²) ) to obtain a modified SCM Cx.; (ii) Prediction: use the SCM C.. and the posterior over the exogenous variables Pujx=x compute the counterfactual distribution Px |X=x
counterfactual is claim about what would have happened if some
Without loss of generality, in this work, we focus on the counterfactual identification (C3) of causal mechanisms fi of multi-dimensional (d > 1) variables X; E X, given only observational data (C1).
3.2 Optimal Transport: Static & Dynamic
Optimal Transport (OT) (Villani et al., 2008; Santambrogio, 2015; Peyré et al., 2019) is a suite of techniques for learning an optimal map between measures that minimises transport cost.
Monge and Kantorovich Formulation. Let P(2)) be the set of probability measures on 5 C R The Monge formulation of OT seeks a map T : Ω - Z between two distributions μr, P(22) that
pushes j forward to v (i.e., Tap- = v/) while minimising transport cost function e: Ω x Ω + I::
Typ=v / c(æ, T(x))
W Wo(pvv v) = inf
e.g. c(z, = [
(2)
However, a transport map T may fail to exist, for instance when is discrete and v is continuous. The Kantorovich formulation of OT relaxes the Monge problem and seeks optimal coupling t E [T4, v), where I(u,,)) v) C P() x 22) denotes the set of couplings with marginals μ and 2²:
[
π arg min rEIl(,.))
c(z, v) dT(x,)).
An optimal IT" always exists, and when a Monge map T' also exists, both formulations coincide. Brenier (1991)'s theorem states that, under fairly general conditions, there is a unique and monotone optimal transport map T = Vá. where ò Rd + R is a convex function. In this work, we show that this result has far-reaching implications for the counterfactual identification of causal mechanisms.
Haaaoou-Brenier Formulation. Benamou and Brenier (2000) showed that the above formulation, known as static OT with quadratic cost, can be equivalently expressed using a dynamic formulation.
In simplified terms, one seeks time-dependent velocity field v Rd [0, 1] + ıd for t E [0, that transports Po to P1 along a flow defined by an Ordinary Differential Equation (ODE):
{ [[]o(XX,, [)13] dt dXx, = v(X,, t) dt, Xo ~ po, (X1)aPo = P1
0* ="=arg arg min VEV
(4)
where Vss is the set of admissible velocity fields, i.e., fields that are measurable, ensure well-posed flow with a unique ODE solution, and yield finite kinetic energy while transporting Po to P1.
Theee are other types of interventions, such as replacing f, by some new mechanism F.. That is. for any measurable sets A, B C Q. we have that (<A ②) =(A) and m(2 B) vBB).
3

---

3.3 Neural ODEs and Flow Matching
Continuous Normalizing Flows (CNFs) (Chen et al., 2018) seek to learn mapping from a simple base distribution Xo ~- Po to the data distribution X1 ~ Pdata using an ODE, whose imm-deepeneett
vector field IRd x [o, 1] - IRd is parameterised by neural network with parameters A:
Vt E [0, 1],
dX,- = vy (X₁; 0) dt,
Pt: = (Xt):Po.
(5)
Flow Matching (FM) (Lipman et al., 2023; Liu et al., 2023; Albergo and Vanden- Ganden-Fijnden, 2023)
trains CNFs simulation-free by regressing a parameterised vector field onto a known target field:
min /. Ex, -Pa Ximpi [X:: (X, |x)]] de, = (1- l)xo + LX1,
(6)
where the target vector field defined be dxx; = oz (X; X1) = (X1 Xo) dt. The simplicity of ot is a result of choosing a simple linear interpolation path (McCann, 1997) between Po and pl. OT has recently been used augment neural ODEs by straightening the sample paths of flows, and speeding up simulations at inference time (Liu et al., 2023; Pooladian et al., 2023; Tong et al., 2024).
4 Counterfactual Identifiability: Theoretical Analysis
We begin by defining the scope of our theoretical analysis in the context of existing counterfactual identifiability results. Detailed proofs for all our theoretical results are provided in Appendix B.
4.1 Problem Statement
Without loss of generality, we focus on the counterfactual identification of the 7th causal mechanism f for multi-dimennional variable X within an SCM C. where dim( = dim(U) = d > 1, given only observational data. Further, we are interested in the Markovian confounding, that is, U is independent of X's parents U IL PAx x as detailed in Definition 4.1 below.
where there is unobserved
Definition 4.1 (Markovian SCM). An SCM is said to be Markovian if Ug II Uj whenever * j. In other words, the exogenous variables in the model are statistically independent of each other, and their joint distribution factorises Pu(u) II. Pe, (u:), where u realisation of U. Thus, Markovian SCMs induce a unique, factored joint observational distribution: (x(x) = -I PX.IPA, (a: pa.).
Counterfactual identifiability (L3, cf. Definition 4.2) results using only observational (C1) data exist for the Markovian case when d = by using a strict monotonicity constraint on f. However, results for the d >> 1 case are still missing, as generalising the monotonicity condition with munuti-dimensional (X, U) variables is non-trivial (Nast-Esfahany al., 2023), and remains unexplored in counterfactual identifiability literature. In the presence of unobserved confounding (i.e. non-Markovian SCMs), Nasr- Esfahany al. (2023) established counterfactual identifiability from observational (L₂) data for: (i) d == 1 given a set of Instrumental Variables (IVs)), I (Imbens and Angrist, 1994); and (ii) d 2 1 given set of variables Z that satisfy the Backdoor Criterion (BC)* (Pearl, 2009), w.r.t. the pair (X, PAx). Definition 4.2 (Counterfactual Identifiability (Pearl, 2009)). A counterfactual query Q is identifiable if for any pair of SCMs ((1) and ((2) we have that Q(((2)) = Q(C(2)), whenever X pel
Nasr- Estahanyy et al. (2023)'s counterfactual identifiability results rely on the bijectivity and mono- tonicity of the causal mechanism f. except for the BC case where bijectivity alone is shown to be sufficient for 2 given that a variability assumption holds. This assumption requires that Z influ- ence U strongly enough for Puiz and its gradient to vary across different values of Z. In Markovian settings, the independence of U and PAx x alone is not sufficient for counterfactual identifiability of f. In addition, bijectivity constraints on ʃ are also insufficient when d > 1, due well-known rotational symmetries of the prior on U (Hyvärinen and Pajunen, 1999; Hyvärinen et al., 2024).
In the following, we characterise the set of constraints and assumptions that permit counterfactual identifiability for multi-dimensional variables in Markovian SCMs from observational data. For this, we connect ideas from (dynamic) optimal transport (Benamou and Brenier, 2000; Villani et al., 2008; Santambrogio, 2015) and graphical causality (Pearl, 2009) to build monotone mechanisms for > 1.
PIVs must be I II and only influence X through PAx.
'If Z blocks all backdoor paths from PAx X; deseenncyyyyyccy of Z is only required for adjustment.

---

4.2 Dynamic Optimal Transport for Counterfactual Identification
counterfactual query pp- {X z, PA pa} can be answered by the deterministic map":
f R R" Vz, pa. pa" : x f(pa". u). where u == f
with pa" denoting the counterfactual parents of æ; and ʃ is invertible w.r.t. the exogenous noise u. Therefore, counterfactual transport map that implicitly abducts the exogenous noise exists:
T : 22th RR R. Ya,, pa, pa" (pa",, pa,, x) f(pa", (pa,
(8)
directly pushing the observational distribution forward to the counterfactual T" TTPPXIPAA X|PAppa
Definition 4.3 (Monotone Operator). A mapping f I Rd -
Rd is monotone in if:
(f u1) [pa, 122), u1 u2) 2 O, Vun, U2 R pa R
Proposition 4.4 (Monotone Counterfactual Transport Map). If mechanism f(pa, u) is monotone in u (Def. 4.3), then the respective counterfactual transport map T (pa", pa, is monotone x. Remark 4.5. In Appendix A we provide gentle motivating example (d 1)-dimension of what monotonicity of 7 entails and why it is important for consistent counterfactual inferences. The subsequent proof provides a generalisation to vectors (i.e., 1) in the sense of monotone operators.
Monotonicity of the counterfactual transport map T w.r.t. guarantees that, for a given intervention on pa. the rank order of factual outcomes is preserved in the counterfactuals. This is important (e.g. for fairness) as it prevents rank inversions across individuals under a given intervention. However, learning the counterfactual transport map directly from data is challenging. since we (almost) never have access to paired or unpaired samples from both the observational and counterfactual distributions.
As we will show, the dynamic optimal transport specialisation of T generalises quantile and rank functions for I, ensuring multivariate rank preservation and consistent counterfactual inferences.
Lemma 4.6 (Unique and Monotone Dynamic OT Mechanism). Let dim (X) dim(U) d > and consider Markovian setting U TL PA. Assume Po and P PA absolutely continuous w.r.t. the Lebesgue measure with strictly positive and bounded densities on bounded, open, convex domains. Let T: IR x R + RRd be the time-l map of dynamic optimal transport flow:
(Te:te 10, [],
dT: = we (u; dt,
T.Po XIPA
P
(10)
Then, there exists a convex function & Rdd - R& such that: T(u; = Vudlu; pa), for Pu(u)-aee u, where is monotone, bijective a.e., and uniquely determined by the pair (Pu, pe (IPA).
Remark 4.7. We draw on Brenier (1991)'s theorem to show that the causal mechanism monotone. However, Brenier's map does not guarantee bijectivity by default. For this, we use Caffarelli (1992)'s standard regularity results to show the OT map is locally diffeomorphic. condition can be satisfied by smoothing empirical target distributions by adding mild continuous noise (e.g. uniform or Gaussian) ensure the resulting distribution admits a smooth density.
practice, the bijectivity
Since dynamic OT maps invertible when the velocity field v, is well-defined (e.g. Lipschitz continuous in space), we can recover the counterfactual transport map by composition (cf. Eq. (8)). By Proposition 4.4, if T(u T(u; pa) is monotone in then T (ppa',,,, pa, 2) is monotone in In the following section, we discuss the indeterminacy induced by the choice of prior on the exogenous noise distribution Pt and characterise the necessary conditions for strict monotonicity of T in x.
4.2.1 Exogenous Prior Indeterminacy
Since the true prior distribution over the exogenous variables Py is typically unknown, one must either rely on domain expertise define prior, choose one for mathematical convenience, or try to learn it directly from data. This choice of prior induces an indeterminacy in counterfactual outcomes.
Next, we show that regardless of the prior Pe choose, under standard regularity conditions, there exists a unique and optimal function g* R >> Id that maps to any other Py we would have chosen.
'ee omit the subscripts for simplicity, ie., PA PAx and similarly for their realisations: pa pa.-
5

---

Definition 4.8 (Transport L3- Cy-Equivalence). Let T0,T(2 R x R - Rd be transport maps
with exogenous priors P p( on IRd, and parent domain R We say they are counterfactually equivalent ~Lg. if and only if there exists bijection g IR + Rd such that 9: p(2) p that is:
T()) ~c T = Ag T(t)) pa)= T (w); pa) for PU -ac. Vpa. (11)
Proposition 4.9 the same counterfactuals and only if they are L3 L---quvialentt in the sense of Definition 4.8:
et al. (2023)). Transport maps T(1, T(2), IRd x R + R produce
Vx, pa. pa" x² T (pa", pa,,)) pa, = (2( (pa", pa,.) > T ~c T
(12)
Remark 4.10. Given only observational data (C1), seek transport map T that is counterfactually equivalent (~e,)) to the true transport map T underlying the data geltaaennrating process.
Lemma 4.11 (Existence of the Prior Transition Map). Let P Pp()) be probability measures Th with finite second moments, both absolutely continuous w.r.t. the Lebesgue measure. Then, there exists a transport map Rdd IRd that is unique p!1) -a.e., monotone, and a.e. bijective.
We are now equipped to prove our main counterfactual identifiability result in the sense of Def. 4.8. Theorem 4.12 (Counterfactual Identifiability in Markovian SCMs). Let dim (X) dim (U) d > 1, and U IL PA. Assume Pu is the continuous uniform measure on [0, 1]2. Let T Rd x R - IR be the time-l dynamic OT map described in Lemma 4.6, which pushes Py forward to P Then, the induced time-I counterfactual dynamic OT map T R22k RR
Red is strictly monotone in x::
'T'"""" p1), (apa", - -T" (paa,, pa, 22). - x2) > 0, Væ,, 22 E TRd pa, pa' E TRA.
(13)
Remark 4.13. This shows a dynamic OT flow mechanism (Lemma 4.6) yields counterfactual optimal transport map I that strictly monotone in æ (i.e., monotonic in the vector sense), thereby extending the classical monotone-quantil notion to multivariate (d >> 1) settings and, under stated Caffarelli (1992) regularity conditions, guaranteeing ~L, identifiability from observational data (L1).
4.3 Non-Markovian Counterfactual Identifiability
4
U
2 - A
A
7
X
Figure 1: Four canonical causal graphs. From left to right: (i) Markovian, (ii) Instrumental Variable, (iii) Backdoor, and (iv) Frontdoor. Dashed bidirected ares denote unobserved confounding.
To ensure wide applicability, we extend counterfactual identification from observational data (C1) to Markovian multivariate
1) settings under the following common criteria (Figure 1):
(i) Instrumental Variable (IV), (ii) Backdoor Criterion (BC), and (iii) Frontdoor Criterion (FC).
Appendix B.3 contains the proofs. Nasr-Esfahany et al. (2023) provided results for only when
BC applies. We expand the IV results to d thanks to the mechanism proposed in Lemma 4.6
being monotone, bijective and uniquely determinable under regularity conditions (Caffarelli, 1992). We also prove ~L, identifiability for FC under similar conditions to BC (Nasr-Esfahany et al., 2023).
5
Counterfactual Transport Maps via Flow Matching
Although bijectivity and variability are sufficient for counterfactual identification under BC settings (see Appendix B.3), Markovian/IV settings require monotone counterfactual OT map to obtain equivalent guarantees (Theorem 4. 1). For the former, we prescribe standard continuous-time flows trained via flow matching, and for the latter, we approximate dynamic OT with these flows.
Counterfactual Inference. Using cedduction-action-preeee for abduction T [z: pa). then forwards time for prediction x" T(u; pa") ). as follows:
Computing counterfactuals with a continuous-time flow model is simple. (Pearl, 2009), we first solve the associated ODE backwards in time
f 0r (21; pa, dt,
x [ ve pa", 0)dd dt,
(14)

---

where pa" is obtained from an upstream intervention (i.e. action) on one or more nodes in the associated SCM. The Picard-Lindelóo theorem states that both initial value problems admit unique solutions if or (2:5 pa, 0) is Lipschitz continuous, a condition that is readily satisfiable in practice.
Markovian OT Coupling. Although theory for Markovian and IV settings is not tied to any particular OT approximation, it remains nontrivial to operationalise at scale with current techniques. To validate our theoretical claims, we offer a flow matching approach with a bespoke Markovian Batch-OT coupling that resolves counterfactual consistency issues in the standard formulation. Batch- OT flow matching (Pooladian et al., 2023; Tong et al., 2024) solves the OT problem on batches of source and target samples to form better pairs, using GPU-enabled solvers (Flamary et al.. 2021), I). This asymptotically recovers the global OT map, but requires large batches in practice, especially in high dimensions (Klein et al., 2025). Notably, Mousavi-Hossee et al. (2025) amortise OT pairing via semidiscrete couplings, reporting improved sample quality with substantially lower pairing cost.
As prescribed by Lemma 4.6, the correct way to formulate OT for Markovian structures must ensure that the marginal Pt remains independent of PA upon coupling |PPAA A naive application of Batch-OT flow matching to a Markovian setting implicitly entangles the parent variables PA and the exogenous noise U, which violates the independence requirement U 11 PA. This happens because upon sampling a random batch of exogenous noise [u(²)/) ~ Pv. and observations ((x()) )) ~ PAPPA solving the OT problem then reassigns each (0() to a pair (e(0),a@() but implicitly couples ( with pa( by association with 26)) This invalidates the abduction process in Markovian settings as U is no longer strictly exogenous, and the counterfactuals will be incorrect.
To fix this issue, we can instead sample batched data from the conditional distribution:
(1()) ) Pu,
pa ~ PEA
(, P PX|PA=pa
(15)
then solve the Batch-OT problem for each fixed value of pa in turn. Simply put, the parent values must match across each element in the batch to satisfy the Markovian constraint that U IL PA.
The objective we optimise can be understood as solving a family of OT flow couplings:
min J/-
E,,22 (1pa)
[12 pa, 0) uy x. pa)
d.P. (pa),
(16)
where the flow is x4 = (1 t)u + tx. Here x (m)) denotes the implicit conditional joint distribution:
1 m m
m pa () 8u)82(1 (x),,
π (m) (u, x pa)
(17)
induced by the OT coupling Tpa E R.mh x RT, for m ii.d. samples drawn according to Eq. (15).
Experiments
We conduct two sets of experiments: (i) constructed scenario where the ground-truth counterfactuals are known and we can, in principle, verify our identifiability claims; (ii) a real-world medical imaging dataset widely used for counterfactual inference. When counterfactual ground truth is not available, we follow Monteiro et al. (2023); Ribeiro et al. (2023) and measure the counterfactual soundness axioms of composition, effectiveness and reversibility. While useful, these metrics alone do not imply identification and should not be construed as evidence of causal validity. For details, see Appendix D.
6.1 Counterfactual Ellipse Generation
For this study, we adapt the counterfactual ellipse generation setup by Easr-Esfahany et al. (2023), Let U E R² be the semi-major and -minor parameters of an ellipse, PA E (0, 2T) be an angle specifying a single point on the ellipse, and X E R² be its cartesian coordinates. The data-generating process is: Pz, u ~ Puizma pa PPAJZ=as and x == (pa, u). By construction, U IL PA Z. and Z satisfies the backdoor criterion (BC) w.r.t PA + X. To induce a Markovian setting, we randomise PA. yielding U 1. PA. The goal is learn a map T (0, 2T)
x 2 + R22 that can infer the set of counterfactual points X = (x* = 7 x:) pa E (0. 2m) which draws the entire ellipse each observed (pa, belongs to. Importantly, since we always know the ground truth X", we can evaluate counterfactual estimates exactly using the mean average percentage error (μAPE).
7

---

Observation (X =x
True CF
Paed CF (NAIVE)
Peed CF (Ous)
1000
7.5
PARS
Marne D.T Paret
PArS Pars
Pare Part
Mara
5.0
2.5
25
5.0
7.5
7.5
EBM (Kau: (.319)
FLOW (JaI. 0.312)
OT-EB (.001)
UTT-FLOWWn (jeac 0.001)
0
-2
Obeervation
CF path
Reversel path (20x)
Figure 2: Counterfactual ellipse generation. (Top) Comparing our OT coupling flow to the naive approach (Section 5). (Bottom) OT maps exhibit greater counterfactual reversibility. counterfactual pa Tpm (z), is reversed by I, pm Tpa- (a*), ) and perfect reversal squares the circle.
BBM FLOw
Table 1: Counterfactual ellipse evaluation. As predicted by our theory, monotone flows are consistent in the Markovian setting. and bijectivity alone is sufficient under the front-door criterion.
6.0
OT ETMM O'T-Rrow
RAPE
2.0
Baselines (Nast Esfahuny al., 2023)
SCHEME
CONDDA
CONDPA,
Markovian
BC2
MARE: (5)
6607
6582
607
1.0
(Ours)
NFE
EBM
FLOW
OT-EBB
OT-FLOW
0.5
2
19.74
19.9::30
1.64: 19
1.601.00
Front door
10
5.79.1.11
5.8011.08
1.42.1.23
1.37.1.14
0.0
50
1.79:t.08
1.672.10
1.42::.22
1,35::.15 18
Stop (x (xP)
2
34.0+ .0
33.6 223
1,21
1.06:: 02
Figure 3: (a) CF error; (b) curl of the vector field during training.
Markovian
10
9.41# 24
9.20: 166
0.95-t.02
0.78:t.01
50
2.32.1.00
2.3011.02
0.93.1.02
0.76.1.01
We build four flow variants: (i) an energy based model flow (EBM), which is curl-free by design; (ii) a continuous- time flow (Lipman et al., 2023) (FLOw)); (iii) an EBM using the family of OT couplings described in Section 5 (OT-EBM); and (iv) same as (iii) but using FLOW (denoted OT-FLoW). To show our OT coupling is necessary for correct counterfactual inferences, we compare with the "naive" (NAIVE) version of a Batch-OT flow, which violates Markovianity (cf. Section 5). For all other details regarding architectures, datasets and experimental setup please refer to Appendix E.
Nasr-Esfahany et al. (2023) reports that their spline flow-based model failed in the Markovian case, with /AP It only succeeded using the BCz scheme (Back-door) with a μAPe of 1% (reproduced at .98 3).. Learning P(X PA) P(X PA, Z) also failed when 2 is a confounder (cf. CONDPA. CONDPA,2) 2). Table shows our flows produce near-exact ground-truth counterfactuals, and using just two function evaluations with OT. Counterfactual reversibility (Monteiro al., 2023) also improved by straighter paths (Figure 2). Further, we experiment with a non-Markovian setting with non-lineat unobserved confounding where the front-door criterion applies. The results validate our theory in that bijectivity alone is sufficient for consistent counterfactual inference in this case.
6.2 Case Study: Chest X-ray Imaging Counterfactuals
To extend our study to high-dimensional settings, we conduct experiments MIMIC-CXR (Johnson al., 2019), widely used dataset for counterfactual inference. Our assumed causal graph follows the baselines (Ribeiro et al., 2023; Xia et al., 2024), and includes SEX (S), RACE (R), AGE (A) and DISEASE (D) variables, where A + D, and PAx = (S, R, A, D} cause the X-ray image X. To parameterise our flow models, we use streamlined version of Dhariwal and Nichol 222 UNet
'Equivalent to Sanchez and Tsaftaris (2021) if the source distribution Gaussian (Gao et al., 2025).

---

Table 2: Counterfactual effectiveness on MIMIC Chest X-ray (192x192). Aucu denotes the absolute
difference in ROCAUC of counterfactuals relative to the observed data baseline. For each variable, our results (blue shade) appear on the right, and baseline results (Ribeiro al., 2023) are on the left. Note the large improvements for e.g. RACE. For more comparisons and ablations, see Appendix F.
SEX (S)
RACE (R)
AGE (A)
DISEASE (D)
BASELINE
AUC (%) f
AUC (%)
MAE (yr)
AUC (%) f
Observed data
99.63
95.34
6.197
94.41
INTERVENTION
AUC (%) 4
| AuC| (%) 4
MAA (yr) 4
|AuUCl (%)
do(S s)
0.370
0.173
11.44
0.583
0.288
0.333
2.490
0.023
do r)
0.070
0.180
8,640
0.050
0.144
0.3941.110
7.010
0.310
do A do D d) (rand)
0.070
0.187
14.64
1.197
0,446
0.8361.08
2.810
0.347..00
0.070
0.067
16.04
0.627
0.371
0.435
3.790
2.280
0.170
0.150=
12.54
0.730
0.300
0.510
0.590
0.640
(dmand Siy malk widie (houtly city sai vibig (doasrd Shy itals wwweeeeeewee........... Ciyy desul kidtichhacty 30s, culs wiined Cweatee D0, mades saite Dieittee Iy male tlacd) Ewaatte Doy with
fe Nakk
(bwwod)
Counterfactual
maki
Reversal
Figure 4: Qualitative counterfactual inference results MIMIC Chest X-ray. We observe faithful, reversible interventions without requiring counterfactual fine-tuning, or classifier(-free) guidance.
architecture (see Appendix F for details). Table 2 and Figure report our main results. We observe substantial improvements over baselines (Ribeiro et al., 2023; Xia et al., 2024) using our flows, across all three counterfactual soundness axioms, and without requiring any costly counterfactual fine-tuning or classifier -feee) guidance. That said, this alone does not imply causal validity. In Appendix F, we report additional comparisons and ablations. We observe performance trade-offs: for instance, OT- FLOw (which assumes Markovianity) outperforms on race interventions but underperforms FLOW on disease interventions, suggesting non-Markovian interaction effects subpar OT approximation. Notably, our Markovian OT coupling substantially improves over the NAIVE OT flow baseline.
Conclusion
Causal claims credible only insofar their identification is defensible, because observationally equivalent models can imply different answers. This work establishes a foundation for counterfactual identification high-dimensional, multivariate outcomes (c.g. images) from observational data. We clarify the conditions under which flow matching yields identified counterfactuals for common causal structures, stating all assumptions and constraints explicitly to invite scrutiny and relaxation. Our results address causal validity concerns in prior work and extend to non-Markovian settings under standard criteria. For Markovian and instrumental variable settings, we characterise a continuous flow that is provably unique and monotone, yielding rank- aankprreerrving dynamic OT map and consistent counterfactuals, Using continuous flows, we demonstrate near-exact ground-truth counterfactuals in controlled scenarios and improved axiomatic counterfactual soundness real images relative to prior methods. Scaling to large problems remains a challenge. We conclude by urging practitioners to assess the causal validity of their estimates when using generative models for causal inference.

---

Acknowledgments and Disclosure of Funding
We thank Aapo Hyvärinen and Charles Jones for helpful insights and discussions. F.R. and B.G. acknowledge the support of the UKRI AI programme, and the EPSRC, for CHAI-EPSRC Causality in Healthcare AI Hub (grant no. EP/Y028856/I). B.G. received support from the Royal Academy of Engineering as part of his Kheiron/RAEng Research Chair.
References
Michael Samuel Albergo and Eric Vanden-Eijnden. Building normalizing flows with stochastic interpolants. In The Eleventh International Conference on Learning Representations, 2023. 2, 26
Joshua Angrist and Guido Imbens. Identification and estimation of local average treatment effects, 1995. 17
Michal Balcerak, Tamaz Amiranashvili, Antonio Terpin, Suprosanna Shit, Lea Bogensperger, Se- bastian Kaltenbach, Petros Koumoutsakos, and Bjoern Menze. Energy matching: Unifying flow matching and energy-based models for generative modeling. arXiv preprint wrxiv:2504.0611 2025. 27
Elias Bareinboim, Juan D. Correa, Duligur Ibeling, and Thomas Ieard. On Pearl's Hierarchy and the Foundations of Causal Inference, pages 507-556. Association for Computing Machinery, 2022. 1, 2, 26
Jean-Davis Benamou and Yann Brenier. A computational fluid mechanics solution to the monge- kantorovich mass transfer problem. Numerische Mathematik, 84(3):375 -333, 2000. 2, 3, 19,
26
Yann Brenier. Polar factorization and monotone rearrangement of vector-valued functions. Com- munications pure and applied mathematics, 14(4):375-417 1991. 2, 3, 5, 18, 19, 21, 25,
26
Kailash Budhathoki. Lenon Minorics, Patrick Blöbaum, and Dominik Janzing. Causal structure-based root cause analysis of outliers. In International conference on machine learning, pages 2357-2369. PMLR, 2022.
Luis A Caffarelli.
The regularity of mappings with convex potential. Journal of the American Mathematical Society, 5(1)999-104, 1992. 2, 6, 19, 21, 22, 26
Guillaume Carlier, Victor Chernozhukov, and Alfred Galichon. Vector quantile regression: an optimal
transport approach. The Annals of Statistics, pages 1165-1192, 2016. 23, 26
Daniel Castro, Ian Walker, and Ben Glocker. Causality matters in medical imaging. Nature Communications, 111)::1-0,, 2020. 1
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud. Neural ordinary differential equations. Advances in neural information processing systems, 31, 2018. 2, 4, 22, 26
Yikang Chen, Lili Tian, et al. Exogenous matching: Learning good proposals for tractable counter- factual estimation. Advances in Neural Information Processing Systems, 37:132113-132363, 2024. 2
Juan Correa, Sanghack Lee, and Elias Bareinboim. Nested counterfactual identification from arbitrary surrogate experiments. Advances in Neural Information Processing Systems, 34:6856-6867, 2021.
2
Fabio Sousa Ribeiro and Ben Glocker. Demystifying variational diffusion models, Foundations and Trends in Computer Graphics and Vision, 1762):76-170, 2025. 27
10

---

Prafulla Dhariwal and Alexander Nichol. Diffusion models beat gans image synthesis. Advances in neural information processing systems, 34:8780-877 87-844 2021. 8, 35
Rémi Flamary, Nicolas Courty, Alexandre Gramfort, Mokhtar Z Alaya, Aurélie Boisbunon, Stanislas Chambon, Lactitia Chapel, Adrien Corenflos, Kilian Fatras, Nemo Fournier, et al. Pot: Python optimal transport. Journal of Machine Learning Research, 22(78):1-55 2021.
David Galles and Judea Pearl. An axiomatic characterization of causal counterfactuals. Foundations 10.102344110990222558444 28
of Science, 3:151-182, 1998. ISSN 1572-8471. doi: 10.
Ruiqi Gao, Emiel Hoogeboom, Jonathan Heek, Valentin De Bortoli, Kevin Patrick Murphy, and Tim Salimans, Diffusion models and gaussian flow matching: Two sides the same coin. The Fourth Blogpost Track at ICLR 2025, 2025. URL https:tt ptttps://operrvvieet 8, 31,33 33
Tomas Geffner, Javier Antoran, Adam Foster, Wenbo Gong, Chao Ma, Emre Kiciman, Amit Sharma, Angus Lamb, Martin Kukla, Nick Pawlowski, et al. Deep end-to-end causal inference. Transactions on Machine Learning Research, 2024. ISSN 2835-8856.
Luigi Gresele, Julius Von Kügelgen, Jonas Kübler, Elke Kirschbaum, Bernhard Schölkopf, and Dominik Janzing. Causal inference through the structural causal marginal problem. In International conference on machine learning, pages 7793-7824. PMLR, 2022.
Joseph Y. Halpern. Axiomatizing causal reasoning. In Proceedings of the Fourteenth Conference on Uncertainty in Artificial Intelligence, UAI'98, pages 202-210, San Francisco, CA, USA, 1998. Morgan Kaufmann Publishers Inc. ISBN 155860555X. 28
Jason Hartford, Greg Lewis, Kevin Leyton- and Matt Taddy. Deep iv: A flexible approach for counterfactual prediction. In International Conference on Machine Learning, pages 1414-1423. PMLR, 2017. 2
Jonathan Ho and Tim Salimans. Classifier-free diffusion guidance. arXiv preprint arXiv:2207.12598. 2022. 34
Patrik Hoyer, Dominík Janzing, Joris Mooij, Jonas Peters, and Bernhard Schölkopf. Nonlinear causal discovery with additive noise models. Advances in neural information processing systems, 21, 2008.
Aapo Hyvärinen and Petteri Pajunen. Nonlinear independent component analysis: Existence and uniqueness results. Neural networks, 12(3):429-439) 1999. 2, 4
Aapo Hyvarinen, Hiroaki Sasaki, and Richard Turner. Nonlinear ica using auxiliary variables and generalized contrastive learning. In The 22nd international conference artificial intelligence and statistics, pages 859-868. PMLR, 2019. 2
Aapo Hyvärinen, Ilyes Khemakhem, and Ricardo Monti. Identifiability of latent-variable and structural couctural-equation models: from linear to nonlinear. Annals of the Institute of Statistical Mathe- matics, 76( 1-33- 2024. 1, 4, 26, 27
Guido W Imbens and Joshua D Angrist. Identification and estimation of local average treatment effects. Econometrica, 62(2):467-475, 1994. ISSN 00129682, 14680262. 4
Alexander Immer, Christoph Schultheiss, Julia E Vogt, Bernhard Schölkopf, Peter Bühlmann, and Alexander Marx. On the identifiability and estimation of causal location-scale noise models. In International Conference Machine Learning, pages 14316-14332. PMLR, 2023.
Adrián Javaloy, Pablo Sánchez-Marttn, and Isabel Valera. Causal normalizing flows: from theory to practice. Advances in Neural Information Processing Systems, 36:58833-58864 2023.
Alistair EW Johnson, Tom Pollard, Seth Berkowitz, Nathaniel R Greenbaum, Matthew P Lungren, Chih-ying Deng, Roger G Mark, and Steven Horng. Mimic-cxr, a de-identified publicly available

---

database of chest radiographs with free-text reports. Scientific data, 6(1):1-8, 2019. 8
Amir-Hossein Karimi, Julius Von Kügelgen, Bernhard Schölkopf, and Isabel Valera. Algorith- mic recourse under imperfect causal knowledge: probabilistic approach. Advances in neural information processing systems, 33:265-277, 2020.
Ilyes Khemakhem, Ricardo Monti, Robert Leech, and Aapo Hyvarinen. Causal autoregressive flows. In International conference on artificial intelligence and statistics, pages 3520-3528. PMLR, 2021.
2,
Michal Klein, Alireza Mousavi-Hosseini, Stephen Zhang, and Marco Cuturi. On fitting flow models with large sinkhorn couplings. arXiv preprint arXiv:2506.05522 05526, 2025. 7, 35, 36, 38
Aneesh Komanduri, Chen Zhao, Feng Chen, and Xintao Wu. Causal diffusion autoencoders: Toward counterfactual generation via diffusion probabilistic models. In ECAI 2024, pages 2516-2523 IOS Press, 2024.
Amar Kumar, Anita Kriz, Mohammad Havaei, and Tal Arbel. PRISM: High-resolution & precise counterfactual medical image generation using language-guided stable diffusion. In Medical Imaging with Deep Learning, 2025. 2
Matt J Kusner, Joshua Loftus, Chris Russell, and Ricardo Silva. Counterfactual fairness. Advances in neural information processing systems, 30, 2017. 1
Yaron Lipman, Ricky TQ Chen, Heli Ben-Hamu, Maximilian Nickel, and Matthew Le. Flow matching for generative modeling. In The Eleventh International Conference on Learning Representations, 2023. 2, 4, 8, 33
Xingchao Liu et al. Flow straight and fast: Learning to generate and transfer data with rectified flow. In The Eleventh International Conference on Learning Representations, 2023. 4
Francesco Locatello, Stefan Bauer, Mario Lucic, Gunnar Raetsch, Sylvain Gelly, Bernhard Schölkopf, and Olivier Bachem. Challenging common assumptions in the unsupervised learning of disentan- gled representations. In international conference machine learning, pages 41 4114-4124. PMLR, 2019. 2
Chaochao Lu, Biwei Huang, Ke Wang. José Miguel Hemúández-Lobato, Kun Zhang, and Bernhard Schölkopf. Sample-efficient reinforcement learning via counterfactual-based data augmentation. arXiv preprini arXivv2012.09092, 2020. 2
Charles F. Manski. Monotone treatment response. Econometrica, 65(6):3311-133 1997. ISSN 00129682, 14680262. 17
Robert.] J McCann, A convexity principle for interacting gases. Advances in mathematics, 128(1):
153-179, 1997. 4, 26
Raghav Mehta, Fabio De Sousa Ribeiro, Tian Xia, Mélanie Roschewitz, Ainkaran Santhirasekaram, Dominic C Marshall, and Ben Glocker. Cf-seg: Counterfactuals meet segmentation. In Inter- national Conference on Medical Image Computing and Computer-Assisted Intervention, pages 117-127. Springer, 2025.
Valentyn Melnychuk, Dennis Frauen, and Stefan Feuerriegel. Partial counterfactual identification of continuous outcomes with a curvature sensitivity model. Thirty-seventh Conference on Neural Information Processing Systems, 2023. 2
Miguel Monteiro, Fabio De Sousa Ribeiro, Nick Pawlowski. Daniel C Castro, and Ben Glocker. Measuring axiomatic soundness of counterfactual image models. In The Eleventh International Conference on Learning Representations, 2023. 8, 28, 31,
Alireza Mousavi-Hosseini, Stephen Y Zhang, Michal Klein, and Marco Cuturi. Flow matching with semidiscrete couplings. arXiv preprint (rXiv:2509.255/9 2025. 7, 35, 38
12

---

Arash Wass-Essfaan and Emre Kiciman. Counterfactual non-identifiabiliy of learned structural causal models. arXiv preprint arXiv:2301.09031, 2023. 26
Arash Nasr-Esfahany. Mohammad Alizadeh, and Devavrat Shah. Counterfactual identifiability of bijective causal models. In International Conference on Machine Learning, pages 25733-25754. PMLR, 2023. 2, 4, 6, 7, 8, 20, 23, 24, 25, 26, 27, 29, 34
Yushu Pan and Elias Bareinboim. Counterfactual image editing. In International Conference on Machine Learning, pages 39087-39101. PMLR, 2024. 2
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, gegh-performance deep learning library. Advances in Neural Information Processing Systems, 2019. 29
Nick Pawlowski, Daniel Castro, and Ben Glocker. Deep structural causal models for tractable counterfactual inference. Advances in Neural Information Processing Systems, 33:857-869, 2020.
35
Judea Pearl. Causal diagrams for empirical research. Biometrika, 82(4):669 6888, 1995. 1,2
Judea Pearl. The logic of counterfactuals in causal inference. Journal of American Statistical Association, 95, 2000. Discussion of Causal Inference Without Counterfactuals by A.P. Dawid.
Judea Pearl. Direct and indirect effects. In Proceedings of the Seventeenth Conference on Uncertainty in Artificial Intelligence. UAT'1,, pages 41 4-10, San Francisco, CA, USA, 2001. Morgan Kaufmann Publishers Inc. ISBN 1558608001. 2
Judea Pearl. Causality. Cambridge university press, 2009. 1, 2, 4, 6, 23, 26, 28
Judea Pearl and Elias Bareinboim. External Validity: From Do-Calculus to Transportability Across Populations, page 452-482. Association for Computing Machinery, New York, NY, USA, edition, 2022. ISBN 9781450395861. 26
Judea Pearl and Dana Mackenzie. The book of why: the new science of cause and effect. Basic books, 2018. 2
Fernando Pérez-Gareía, Sam Bond- Taylor,ylor, Pedro Sanchez, Boris van Breugel, Daniel C. Castro, Harshita Sharma, Valentina Salvatelli, Maria T. A. Wetscherek, Hannah Richardson, Matthew P. Lungren, Aditya Nori, Javier Alvarez- Valle, Ozan Oktay, and Maximilian Ilse. Radedit: Stress- testing biomedical vision models via diffusion image editing. In Computer Vision - ECCV 2024, pages 358-376, Cham, 2025. Springer Nature Switzerland. ISBN 978-3-031--325446. 42
Jonas Peters, Joris M Mooij, Dominik Janzing, and Bernhard Schölkopf. Causal discovery with continuous additive noise models. The Journal of Machine Learning Research, 15(1):2009-2033, 2014. 2
Jonas Peters, Dominik Janzing, and Bernhard Schölkopf. Elements of causal inference: foundations and learning algorithms. The MIT Press, 2017.
Gabriel Peyré, Marco Cuturi, et al. Computational optimal transport: With applications data science. Foundations and Trends in Machine Learning, 116-6))3355-07,, 201, 2,
Silviu Pitis, Elliot Creager, Ajay Mandlekar, and Animesh Garg. Mocoda: Model-based counterfac- tual data augmentation. Advances in Neural Information Processing Systems, 35:18 181566 2022. 1
Aram-A Aram-Alexandr Pooladian, Heli Ben-Hamus Carles Domingo-Enrich. Brandon Amos, Yaron Lipman, and Ricky TQ Chen. Multisample flow matching: Straightening flows with minibatch couplings. In International Conference on Machine Learning, pages 28 1800-28812.. PMLR, 2023. 2, 7, 27,
35, 36
13

---

Rajat R Rasal, Avinash Kori, Fabio De Sousa Ribeiro, Tian Xia, and Ben Glocker. Diffusion counterfactual generation with semantic abduction. In Forty-second International Conference on Machine Learning, 2025. 2
Fabio De Sousa Ribeiro. Tian Xia, Miguel Monteiro, Nick Pawlowski, and Ben Glocker, High fidelity image counterfactuals with probabilistic causal models. In International Conference Machine Learning, pages 7390-7425. PMLR, 2023. 2, 7, 8, 9, 28, 35, 36, 37, 39, 40
Mélanie Roschewitz, Fabio Sousa Ribeiro, Tian Xia, Galvin Khara, and Ben Glocker. Counterfac- tual contrastive learning: robust representations via causal image synthesis. In MICCAI Workshop on Data Engineering in Medical Imaging, pages 22-32. Springer, 2024.
Pedro Sanchez and Sotirios A Tsaftaris. Diffusion causal models for counterfactual estimation. In First Conference on Causal Learning and Reasoning, 2021. 2, 8, 33
Pablo Sánchez-Martnn Miriam Rateike, and Isabel Valera. Vaca: Designing variational graph autoencoders for causal queries. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 36, pages 155--168, 2022. 2
Filippo Santambrogio. Optimal transport for applied mathematicianss, volume 87. Springer, 2015. 2, 3,4
Bernhard Schölkopf. Causality for Machine Learning, pages 765-804, Association for Computing Machinery, New York, NY, USA, edition, 2022. ISBN 9781450395861.
Bernhard Schölkopf, Francesco Locatello, Stefan Bauer, Nan Rosemary Ke, Nal Kalchbrenner, Anirudh Goyal, and Yoshua Bengio. Toward causal representation learning. Proceedings of the IEEE, 109(5)x612-334 2021.
Shohei Shimizu, Patrik O Hoyer, Aapo Hyvärinen, Antti Kerminen, and Michael Jordan. A linear non-gaussian acyclic model for causal discovery. Journal of Machine Learning Research, 7(10), 2006. 2
Ilya Shpitser and Judea Pearl. What counterfactuals can be tested. In 23rd Conference on Uncertainty in Artificial Intelligence, UAI 2007, pages 352-359, 2007.
Jin Tian and Judea Pearl. Probabilities of causation: Bounds and identification. Annals of Mathematics and Artificial Intelligence, 28(1):287-313 2000. 2
Alexander Tong, Kilian Fatras, Nikolay Malkin, Guillaume Huguet, Yanlei Zhang, Jarrid Rector- Brooks, Guy Wolf, and Yoshua Bengio. Improving and generalizing flow-based generative models with minibatch optimal transport. Transactions on Machine Learning Research, 2024. ISSN 2835-8856. 2, 4, 7, 25.
Stratis Tsintsis and Manuel Rodriguez. Finding counterfactually optimal action sequences in con- tinuous state spaces. Advances in Neural Information Processing Systems, 32:3220-227,, 2023. 2023
Cédric Villani al. Optimal transport: old and new, volume 338. Springer, 2008. 2, 3, 4, 19
Athanasios Vlontzos, Bernhard Kainz, and Ciarán M Gilligan-Lee. Estimating categorical counter- factuals via deep twin networks. Nature Machine Intelligence, 5(2):159-166, 2023. 2
Yulun Wu, Louis McConnell, and Claudia Iriondo. Counterfactual generative modeling with varia- tional causal inference. In The Thirteenth International Conference on Learning Representations, 2025. 2
Kevin Xia, Yushu Pan, and Elias Bareinboim. Neural causal models for counterfactual identification and estimation. In The Eleventh International Conference on Learning Representations, 2023. 2, 26
14

---

Tian Xia, Mélanie Roschewitz, Fabio De Sousa Ribeiro, Charles Jones, and Ben Glocker. Mitigating attribute amplification counterfactual image generation. In International Conference on Medical Image Computing and Compuler-Assisted Intervention, pages 546-556. Springer, 2024. 8, 9, 35, 36
Tian Xia, Fabio De Sousa Ribeiro, Rajat R Rasal, Avinash Kori, Raghav Mehta, and Ben Glocker. Decoupled classifier-free guidance for counterfactual diffusion models. arXiv preprins arXiv:2506. 14399, 2025. 2
Junzhe Zhang and Elias Bareinboim. Fairness in Deeisiion-maiin----a the causal explanation formula. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 32, 2018.
Kun Zhang and Aapo Hyvärinen. On the identifiability of the post-nonlinear causal model. In Proceedings of the Twenty-Fifth Conference on Uncertainty in Artificial Intelligence, pages 647- 655, 2009. 2
15

---

Appendices
Table of Contents
A Motivating Example: Monotonicity Requirement
16
B Proofs: Counterfactual Identifiability
18
Optimal Transport for Counterfactual Identification
B.2 Exogenous Prior Indeterminacy
20
B.3 Non-Markovian Counterfactual Identifiability
23
B.4 Assumptions: Summary & Plausibility
26
C Energy Based Models: Curl-free Flows
27
Counterfactual Soundness Axioms
28
E
Counterfactual Ellipse: Extra Results
29
E.1 Dataset Details
29
E.2
Architecture and Experimental Setup
29
E.3
Qualitative Results
30
E,4
Ablation Studies: Priors & Assumption Violations
................................
34
F
MIMIC Chest X-ray: Extra Results
35
F.1
Dataset Details
35
F2
Architecture and Experimental Sctup
35
F3
Additional Comparisons and Ablation Studies
36
F4
Qualitative Results
39
A Motivating Example: Monotonicity Requirement
We start with a simple instructive example in one dimension (d -d=1 1). Suppose we have two causal mechanisms (1) and (2) defined as follows:
T (PA, U) PA + U,
(18)
X T U) PA + 1 1-U,
(19)
where PA ~ ((0.5) is Bernoulli distributed, and U ~ ((,0 1) is continuous uniform. One can verify that both (1) and (2))) induce the same conditional observational distribution:
P(X PA 0) u(0, 1),
P(X PA 1) 211
(20)
Without loss of generality, consider the counterfactual query
X xp (PA X 0.8).
(21)
In the abduction step, (1) and (2) infer different exogenous noise values:
n (T( (pa,
pa 0.8 0 -= 0.8, (2) (7( (pa, = pa +
(22)
0+10+1100 0.8 0.2.
(23)
Nonetheless, both mechanisms produce the same counterfactuals:
x (1) = 7 (pa", ,u(")) ppaa + (2)) = 0.8 1.8,
(24)
r (2) = 7 (pa", 12) -pa" +1 n = =
0.2 = 1.8.
(25)
16

---

Even if () # u(2), both (()) and 7(2) can return the same counterfactual outcomes because they assign the same rank u (i.e., cumulative probability or quantile level) to the observed value æ.
Concretely, since 7(1, T define the same observational distribution P(X PA), their conditional CDFs FXIPA R + [0, match, and we can use them to perform abductionn:
u = FXIPA-0(2)
(26)
= P(X < x PA = =0)
(27)
0.8 - 0
0.8,
(28)
b-a a
1-0
where (a, b) = (0, 1) since P(X PA = 0) = =1(0,1).
Then, we compute the counterfactual using the quantile function F=1 the inverse conditional CDF under intervention:
X/pPAppa' (0,1) that is.
EXIPAppi pa (u) = a a u(b a) 11+ -18,
(29)
where (a. b) (1, 2), since we set pa" == and P(X PA ()1
The key here is the strict monotonicity T(.), (()) in U, that is for any (1)) < (2) E IR we have
T (1, <T (.))())
and
(., ( v22)).
(30)
To illustrate the importance of this property, consider third mechanism (3) that induces the same conditional observational distribution:
U if =
X = = (PA,))
U ~ U(011)
(31)
2
if =
However, we can see that (3) produces a different counterfactual for the same query used before:
u (3) (T(3) (pa, r) = x = 0.8,
(32)
x (3) = (r(3 (pa", u (3)) = 2 - v = 1.2.
(33)
One can verify that, unlike for T()) 72(2), the rank is not being preserved under intervention:
ox) = 0.8 > FXIPA-pa- (xPPA-pa-(I*(3)) = 0.2.
(34)
For a consistent counterfactual inference model, we would expect the inferred rank u (i.e., the cumulative probability or quantile level) of an observed outcome x given its parents pa to map to counterfactual outcome x² that lies the same rank the counterfactual distribution:
F F
Fx|PA=pa (Fx PA))a
(35)
The mismatch observed previously occurs because (3) monotonic in U across different values of the parents PA". As result, (3) fails to preserve the rank order of outcomes across different interventions, thereby breaking counterfactual consistency and identifiability.
This pank-preservation condition can be viewed as the continuous analogue of the monotone trealment response (Manski, 1997; Angrist and Imbens, 1995) used in epidemiological studies, which rules out so-called 'defiers': units for whom exposure to discase risk factors would actually reduce disease risk. both cases, monotonicity ensures that counterfactuals respect the domain knowledge- based intuition that exposure should not make disease less likely. As is often the case in causal analysis, the appropriateness of this condition contingent on the particular context in which it is applied.
With the above d = example in mind, in the following sections, we characterise the set of constraints on the map T' which permit (point-wise) counterfactual identification of causal mechanisms for multi- dimensional variables (i.e. d d> 1); a key gap in existing literature. To achieve this, we connect ideas from dynamic optimal transport theory and graphical causality in ways previously unexplored.
Noot that (3) is monotonic in U/ for each fixed value of PA. but not globally across both settings.
17

---

B Proofs: Counterfactual Identifiability
B.1 Optimal Transport for Counterfactual Identification
Proposition B.1 (Monotone Counterfactual Transport Map). If mechanism æ = fpa, u) is monotone in u, then the respective counterfactual transport map T (pa", pa, x) is monotone in x.
Proof. If f(a, u1) < f fpa, u2) for all U1 < U2 then its inverse must be monotonic in given pa: (pa, 222)
121.22
22 f (pa, æ1)
(36)
Applying the same logic given pa". we have that their composition is also monotonic 2:
12 21 x2 = T (pa", pa, T (pa", (2)
(37)
To generalise to multiple dimensions (d 1), assume ʃ R x Rd - IRd is monotone in w:
Vui, (f(pa, u1) f(pa, u2). 11 U2n N 0,
(38)
and let
3 x2 = (pa, u1) - (pa, u2),
"1 - 1/2 =1
(pa, - 1000-s--pppp x2).
(39)
By substituting in the above, we have that:
(/(pa, u1) - [paa, w2), 01 u2) = (I) - x2, / (pa, - -5- (pa' 2 0
(40)
which by symmetry of the dot product, shows f-1 is monotone in Lastly, to show when the counterfactual transport map 7 (ppa", ,a,, x) - /paa / (pa, x)) is also monotone in x (for d >>>> 1), suppose that the causal mechanism is equal to the gradient of convex function:
R R
+E E
pa f(pa, u) Fup(u:pa).
(41)
Recall that gradients of convex functions are monotone:
Yu U1, ug pa Vuo(uu;pa) V pa), u1 U2) 2 0,
(42)
and therefore, the following must hold:
V yx pa, pa" (ppa", pa, x) = Vuo(u; pa paa) ((aaa))
(43)
V2 (T* (pa", pa, x1) T (pa', pa, x2), 21 - x2) 20, 0,
(44)
where (v. w)2 = v 00 (pf x x),, concluding the proof.
Remark B.2. Monotonicity of the counterfactual transport map T w.r.t. x guarantcess that, for a given intervention, the rank order of factual outcomes is preserved in the counterfactuals. This is important for, e.g. fairness, as it prevents rank inversions across individuals under intervention. Euclidean monotonicity of T (for d > 1) is guaranteed if the symmetric part of its Jacobian 0, fu,ff is positive semi-definite, this holds, for example, if the two Jacobian components commute.
Next, we present Brenier's theorem, celebrated result in optimal transport theory that we draw from in our counterfactual identifiability analysis of dynamic OT-based causal mechanisms.
Theorem B.3 (Brenier (1991)). Let μ. v be probability measures on Ω C Rd with finite second moments, and suppose is absolutely continuous w.r.t. Lebesgue measure. Then, for the quadratic cost function c(x, T(x)) = || - T(x) the Monge problem admits a unique optimal transport map T Ω. Q pushing μ forward convex T Vó, and is monotone in the sense that:
i.e., Tel V. Moreover, this optimal map is the gradient of
Væ1.32 (T(x1) T(x2),
0.
(45)
Remark B.4. This result provides powerful geometric insight: under regularity conditions, the most efficient way of transporting mass from a source to a target distribution is by pushing along the gradient of a convex function. The gradient structure ensures that the map does not 'fold', and the monotonicity condition reflects directional consistency in the sense that points which are initially close tend to remain close under the map, thereby preserving order and preventing overlaps.
18

---

Lemma B.5 (Unique and Monotone Dynamic OT Mechanism). Let dim(X) = dim(U) = d 1, and consider the Markovian setting U ll. PA. Assume Pt and Px PA are absolutely continuous w.E.f. the Lebesgue measure with strictly positive and bounded densities on bounded, open, comvex domains. Let T R² x R IRd be the time-l map of dynamic OT flow satisfying:
{Ti e [0,
du; = (TT((u pa)) dt,
T.Po P
(46)
Then, there exists a convex function Φ: - - R. such that: T(u: pa) Vub(u: pa), for Pu (u)-)a.e u, where T is monotone, bijective a.e. and uniquely determined by the pair (Pu, P (PA)
Proof. First, we invoke standard Caffarelli (1992) regularity conditions for OT (see Theorem 12.50 in Villani et al. (2008) for helpful summary). Let 22 A C Rd be bounded, open, convex sets. The source Py and target X IPA are absolutely continuous w.r.t. Lebesgue measure, with densities Pt and P PXIPA bounded above and below by positive constants on 12 and A, respectively. These conditions cover the strictest version of our results; weaker conclusions follow under standard relaxations.
According to Benamou and Brenier (2000)'s ' dynamic formulation of optimal transport, the time-1 map T of a flow (T, L E [0, }]} minimising the kinetic energy
inf (.) (pe
J J [Iu(uupa)]] (e pa) du dpadt
(47)
subject to the continuity equation
Qept + div(prv))
Po PU:
PI PXIPA
(48)
solves the Monge problem:
/ (u: 'd
7 arg min T Po-PA),P==PX
T(u;
(49)
Since Pu < Ld (i.e. absolutely continuous w.r.t. the Lebesgue measure), Brenier (1991)'s theorem (cf. Theorem B.3) applies, and implies that the optimal transport map T is of the form:
T (w; pa) = V Vuo(u;aa)
(50)
where Φ is convex on Ω. Further, this gradient map is well-defined Pv-almost everywhere and is uniquely determined by Py and pe PPI Moreover, since V pushes Pt forward to pe PAIPAA this implies that for almost every u E Q2 it satisfies the Monge-Ampère equation:
(u) P (Vno(u; pa)) det Dio(upa)) pa).
(51)
If is differentiable and strictly convex, then V is strictly monotone and therefore injective:
(Vnó(u(); pa) v.)(((²); pa), (1) u2))
for (1))
+ (2))
(52)
In particular, V Vn pa) (2)( pa) implies (1) 22) and thus is injective almost everywhere. Moreover, since Till- pe PA no Borel set A c A of positive pe PPP -measure can be missed by T and must satisfy Po (T (A)) > O, hence T is onto A up to null sets.
The regularity of & now follows from the theory of elliptic partial differential equations, specifically from Caffarelli (1992)'s 's well-known interior and global regularity results for convex solutions the Monge-Ampère equation. If pu and px PA are bounded above and below by positive constants on
PA
D. A and are C* then by Caffarelli (1992)'s theory the convex potential á E C (2); in particular, Di u; pa) is almost everywhere positive definite and Vo locally a diffeomorphism. Hence, we conclude that T(u pa) = V V. pa). Po-almost everywhere, and that T is monotone, almost everywhere bijective, and uniquely determined by the pair (Pu, PC P
Remark B.6. We draw on Brenier (1991)'s theorem to show that the causal mechanism is monotone. However, Brenier's map does not guarantee bijectivity by default. For this, we must invoke Caffarelli (1992)'s standard regularity results for OT, to ensure the map is locally diffeomorphic. In practice, this condition often satisfied by smoothing empirical target distributions by adding mild continuous noise (c.g. uniform or Gaussian) to ensure the resulting distribution admits a density.
19

---

B.2 Exogenous Prior Indeterminacy
Given only observational data (L1), our goal is to find a transport map T that is counterfactually equivalent (~ec,) to the true transport map "* underlying the data tatt--eenraatig process.
Definition B.7 (Transport Lg-Eqquivalence) Let T(1, T(2) IRd x R* + Rd be transport maps with exogenous priors P pl on RR, and parent domain R We say they are counterfactually equivalent ~Lgi if and only if there exists bijection g: : IRd Rd such that g -P that is:
T()) Le (:22) + 3g : T pa) T (w); for -a.e u, Vpa.
(53)
Proposition B.8 the same counterfactuals if and only if they are eg-eqivvalent in the sense of Definition B.7:
et al. (2023)). Transport maps T(1, (()) R x XR* + R produce
Vx, pa. pa' x² = = (pa", pa, x) = 7 (2) (pa" pa, + T() ~2s T2)
(54)
Proof. We provide a proof here for completeness. An arbitrary counterfactual query Xpa {X- = x, PA pa), where {X x, PA pa} is the observed evidence, answerable by any two equivalent transport maps T(1), T22) IRd x R Rd via aboluction-actioo--rediie
We first prove $
In the abduction step, both maps use the observed evidence to infer the exogenous noise:
(1) (T(1)) (r;; pa),
u (2) (T(2)) (2;; pa),
(55)
and their equivalence (cf. Definition B.7) implies that:
3g ( - g(u
(56)
In the prriin-preiittion step, for any intervention on the parents, we have that:
pa") = pa") thus both maps produce the same counterfactual outcomes:
TT)
(g( pa") = T-(2 (g- pa") = (x(2), pap))aa))))
(57)
Vx, pa. pa" : x²= =T=''''ppa" pa, x)) = T (pa", pa,.). æ.)
(58)
We now prove +
Suppose, for the sake of contradiction, that the function g depends on the parents pa, then:
())) pa : (" pa) = T (g (u pa)::aa)
(59)
Since T(1, T(2) ((2) produce identical counterfactuals, they must agree for any fixed pa":
T paa) = (g (1) pa): pa); pa"),
(60)
where the abduction step is given by:
g (nl)) pa) = (2) (T(2)) -1 pa:
pa).
(61)
Since Eq. (59) also holds when pa = pa", we have that:
T pa") = ((2) ((x()) pa"); pa")
(62)
(2 (1) pa): pa") (2) ((u pa"); pa").
(63)
Since ((2) is injective (i.e. bijective in u (2) by construction). the above implies that:
Va pa" g (2) pa) (2 pa",
(64)
pa g (2 pa) - ('
(65)
This contradicts the assumption that depends on pa. completing the proof.
20

---

Lemma B.9 (Existence of the Prior Transition Map). Let P P(()) be probability measures on I& with finite second moments, both absolutely continuous w.r.t. the Lebesgue measure. Then, there exists a transport map R - Re that is unique P! aae., monotone, and a.e. bijective:
(( inf 91 9:P0 =P
duu aP (u),
c(u.g(u) = 1|u -llllggle =1]=-g(w\]1
(66)
Proof. If PU P have finite second moments and (()) & L by Brenier (1991) (Theorem B.3),
there exists a unique P pf))a.c.. -ac.. monotone, map g = V.. Since Pi < <L', the map is a.c. bijective. Under standard Caffarelli (1992) regularity conditions (bounded densities and convex domains), between supports.
can upgrade a.e. invertibility to a homnomorppisn/diieeomppssss
To fix common reference, we invoke the higher- dimensioall Probability Integral Transform (PIT), known as the Rosenblatt transform. It states that any absolutely continuous probability measure Pt on R& (with d > 1) whose density is strictly positive almost everywhere can be represented as the push-forward of the uniform ([0) 1)4) via a measurable and almost everywhere invertible map.
More concretely, first consider the one one-dimensional case (d = 1). For any random variable X E R.
with cumulative distribution function (CDF) Fx + [0, 1]. the PIT states that:
U == Fx)).
X FU)
(67)
where is uniformly distributed on [0, I]. To extend this to higher dimensions (d > 1), we use recursive component-wise construction known the Rosenblatt transform, explained below.
Let U (4(0 1] and define map T [0. 1] Fd component wise in a triangular manner. Further, write (#1;-.., : ua) for coordinates in [0, 1j and let X = (X1,...,X, Xa) be random vector distributed according to Pu. Now set Ti (2z) Fx (uz).
For each k = ., let Xx be the conditional CDF of X given (X1,... then define:
XA-1),
T Under mild continuity assumptions each conditional CDF (which follow from strict positivity of the density of Pt),, this construction is well-defined for almost every (uI...., ua) E [0, 1].
Ta [u up) = [x
1k)1 (ux). (68)
The resulting map
T(u ua) = (Ti(u,), T2 u2), .... Ta Ta(uu..... ua))
(69)
is measurable and invertible almost everywhere. Since U ~ ~(( 11 one can verify that T(U) has law Pv, or equivalently that (10, 1] Pv. Finally, for any two distinct absolutely continuous
probability measures ((2) and pf2) on Bdd one may similarly construct maps T(1), T(2) [0, 1)
Rd such that:
T 149 - p())
T,0, =P(?)
(70)
With an elerseeeeerywhhree inverse (1(() (i.e. (T(2)) - T2)) id a.c.), it follows that the transport map
g = (TT(2) T
satisfies
9: P 2())
(71)
Thus, both measures p(() and p((2) can be considered forwads of the same uniform measure on [0, 1jd, via different invertible maps, showing that any particular choice of Pu is essentially arbitrary. Without loss of generality, we may therefore fix the exogenous prior as Pu = ([00 1)d). and interpret all subsequent distributions as push-forwards thereof, concluding the proof.
Remark B.10. We use the Rosenblatt transform only to fix common reference for the result; any absolutely continuous law can be written as (10, 1]d) for measurable, i.e.-invertible triangular T. This does not imply that T' is the quadratic-cost OT (Brenier) map. Note that the Rosenblatt transform depends on the chosen coordinate order, meaning that different permutations yield different triangular maps (each still pushes ((( 1j4) to the same target law). This indicates that rank-preserving continuous- time flows without OT are realisable by, for example, assuming a coordinate order and parameterising the flow with component autoregressive model. By contrast, the Brenier OT map is invarint and is determined solely by the two measures and the quadratic cost.
21

---

Theorem B.11 (Counterfactual Identifiability in Markovian SCMs), Ler (im dim(U) = d> 1, and U 11. PA. Assume Py is the continuous uniform measure on [0, 1]. Ler Rd x R -, R be the time-1 dynamic OT map described in Lemma 4.6, which pushes Pu forward to P X Then, the induced time-1 counterfactual dynamic OT map T 1B R +R is strictly monotone in x:
(pa", pa, x1) - T (pa", ,a, T2), 21 - r2) > 0, Væ E IR pa, pa E ER.
(72)
Proof. We present generalised notion of the monotonicity requirement for mechanisms (d > 1), and show counterfactual identifiability can be achieved using the dynamic OT map described in Lemma B.5.
causal
2 and U 1 PA. Let Py be a uniform measure U([0, 1]"). and assume PC IPA is absolutely continuous w.r.. the Lebesgue measure. Let T Rd x R - TRd be the time-1 map a dynamic OT flow {Tr : [0, 1]} described Lemma B.5, that is:
Let dim (X) dim
dT, pa) ex (1; pa)) dt,
T(: P(:PA),,P P PA-
(73)
By Theorem B.3, there exists an optimal map that minimises quadratic transport cost:
T(u; pa) = V.o(u; pa),,
(74)
where o: Rd + R. is a convex function. Importantly, T is well-defined Pu-almost everywhere and is uniquely determined by Py and pe PPA This gradient map generalises the notion of miedimmenional monotonicity, as it is 'vector monotonic' in the sense that:
Vu1, U2 pa (T(u1:pa) - T(u2: pa), d1 (2) 2 O.
(75)
Since Pe P|PA assumed be absolutely continuous, under Caffarelli (1992) regularity conditions (Lemma B.5), the map also bijective. Then, follows that:
T (T (u: pa): pa) u,
[0, 1J
(76)
where the respective maps are realised by solving the associated ODE backwards in time for abduction:
L. ee (T (u: pa) dt,
Ug T (z; pa)
(77)
and forwards in time for prediction of x², for instance, under chosen intervention pa":
T(u:: pa*) uu + L. ve (TT pa")) dt.
(78)
By the Picard-Lindellff theorem, both initial value problems admit unique solutions if the associated velocity field is Lipschitz continuous, which is readily satisfiable in practice (Chen et al., 2018).
Under Caffarelli (1992)'s standard regularity and smoothness conditions (Lemma B.5), the time-1 I
dynamic transport map T arises as the gradient of a strictly convex potential T(u = Vuo(w; pa). and its Jacobian a, T (u; pa) Vig(u, u; pa) is symmetric dositive-definit (SPD).
The Jacobian of the counterfactual transport map
THd (pa", pa. "ppa",pa.u
== T(T (z; pa); pa")
(79)
is thus given by the product of two SPD factors:
0 (pa", pa. x) a. (w: pap) I 0 (x:
(80)
For two distinct points T1, T2 E Ω. define the straight line segment
V1 E [0, 1], where: = + 0. Applying the fundamental theorem of calculus, then the chain rule, we get:
Tt (1 l) + txi, dx, (x) x2) di = dl,
(81)
(ppp",, pa...)) 21) Tppp",, pa, x22) S''' d
dt (' pa, m2)] df
(82)
[ On,T" pa, dt.
(83)
22

---

Taking the inner product with 2 in the form (e, w). == (x:: pa) w yields:
{" (0s. T" 2) dl,
(P'"p,,, pa, - TT' (pa", pa. 02), =
(84)
! 2 10TTT 2d
(85)
L (0,, (0,,T~>>) T 0,T (0,, 2) di > 0,
(86)
because 0,T is SPD and 0,,T~'= + 0 for each E [0, 1] (since 0, SPD). Hence the integral
is strictly positive and the counterfactual transport map T is strictly monotone in as required. D
Remark B.12. An intuition for the above result is provided next. Since the time-1 transport map
T: [0, 1] x I + Rd exists, is uniquely optimal, and (Brenier- monotone, it generalises the notion of a scalar quantile function inverse map 7-1 : Rdd x RR - [0, 1] also exists and shares these properties, it serves as a vector rank function, and can be used to uniquely recover analogously the d case in Appendix A.
vector-valued outcomes (Carlier al.. 2016). Further, since the
The proof above establishes that 7 strictly monotone in æ, and this represents the multivariate (Brenier/OT) notion of monotonicity and is sufficient for rank preservation (Carlier et al., 2016). If one additionally wants coordinate-wise (product-ordes) monotonicity, it holds under extra structure. For example, when Vi pa) has nonnegative aff-diagonal entries, so [:) is coordinate- wise increasing. If, for both pa and pa', , the map T:: pa) is an order-isomorphism for the prod- uet order (hence T (-: pa) is also coordinate- wise increasing), then the composition of the two T" T"pa"",a,
pp,",, x) = T(T- (x; pa); pa") )) is coordinate-wise increasing in œ.
Importantly, for quadratic cost with absolutely continuous source, the OT map is unique a.e. and does not depend on a coordinate ordering; monotonicity in the Brenier/OT sense (i.e., is the gradient of convex potential). This sidesteps the non-uniqueness inherent in coordinate-wist (product) orders, or Rosenblatt transforms, which depend on an arbitrary permutation of coordinates.
B.3 Non-Markovian Counterfactual Identifiability
In the following, extend our theoretical analysis to non-Markovian settings. We provide coun- terfactual identifiability results from observational data under the following standard criteria: (i) Forward Criterion (FC); (ii) Instrumental Variable(s) (IV); (iii) Backdoor Criterion (BC). The proofs for BC (d 2 1) and IV (d = 1) criteria are provided by Nasr-Esfahany et al. (2023). We provide a proof for FC (d 2 1) and extend Nasr-Esfahany et al. (2023)'s IV result d >>> 1.
B.3.1 Frontdoor Criterion Setting
A setting in which counterfactual identifiability under an unobserved confounding U on (A, X) is possible (under stated assumptions) is through the frontdoor criterion, using mediator variable Z (Pearl, 2009). For FC to apply, the following must hold: (i) Z intercepts all directed paths from A to X; (ii) no backdoor path from to Z exists; (iii) blocks all backdoor paths from Lg-equivalence result presented below also applies to the multi-dimensional (d 1) setting.
The
Lemma B.13 (FC Eguiquival,,,e) Consider two FC models T and t with the same observa- tional joint law (A,2,X), where + Z + X and the pair (A, X) is confounded by U, and Ü, respectively. The bijective structural assignments are
T(w; (Tz(a). (x)))) and (fz(a) Tx(a)
(87)
The maps T and I are counterfactually equivalent (~rc)
1. (Frontdoor Criterion):
U LL Z A and t LL LI X AL AI and X 11 A (Ū, Z}.
(88)
2. (Regularity) For every z, the derivatives V. det Jz (+:2) and V det
/et det Jjy (+3z) exist.
23

---

3. (Variability) For every there exist d + 1 points ai...., ad+1 such that
Pula(u an)
V VuPoja(u an)
M (u; .i,,,,..., 0d+11)
(89)
Puja(u a041)
VuPuia(u aa+t)]
is full rank, ie. det M +0. O.
Proof. If Z A is deterministic in both models, then matching the joint law of (A, Z) forces:
2 =AA
2 2=
(2(A) (A) a.s.
(90)
If A is discrete with P(A = a) > 0, or has density strictly positive on an open A and Tz, Zz are continuous on A. then (a)a) = [a)(a) for all If Z A is instead stochastic, with bijective maps w.r.t. noise variables W, W. counterfactual equivalence (Definition B.7) states that:
Tz ~c Tz + 3h, Vw, a Tz(w; a) = f2 -T)(h=(ww)))) (w)
(91)
where h is bijection. By Lemma B.9, h exists, therefore, Tz ~ îz holds.
Next, we prove Tx ~e Tx under similar conditions Nasr-Esfathe Easr-Esfahany al. (2023)'s BC result. First, define link function of both transport maps on the exogenous noise channel:
g(: Tx (xxx((12)) (1:2) (x)+;2),
u g(;; 2).
(92)
which, by the change-of vaaiiallf-variaalls variables formula:
PūJA, (á a. 2) = PUlA,Z (g(u; (i;; 2) a, [a det J (ú: 2
(93)
Under FC (condition 1), we have that U 1L Z A and Ú IL Z A. respectively:
PoJA (0 a) Pila (g( 2) a) det J ( (ū; =) Given regularity (condition 2), differentiate both sides of Eq. (94) w.r.t. 2j:
(94)
0 = (VuPoJA (u a)) 0g(û;; z)
a Oz.
det Jet.J [ + PolA (u a)
det J (ū;
(95)
0z:
and write block form:
8
det J (ú; z)
(pus(ulaa) (VaPola(u a)) a)) []
Uz.
0.
(96)
0g(ū:z)
det J (ūi; 2)
dz:
Stacking the respective equations for a1 .2...... ... ad+1 then yields:
det J (ū;
M(u al,..., ad+1) (d:z)
= 0,
(97)
(g()) dct Jy (;; 2)
U2;
because M must be full rank by variability (condition 3) have that, for all indices i {1,.... d). 0g(û; 22, and 0/0z, det Jy (ū; z) 0. Therefore, g(û; 2) does not depend 2, and can be written as g(ū; 2) = g()). proving Tx and T' differ only by a bijection of the exogenous noise. Finally, combining the above ~c, - equivalencee results, there exist bijections and h such that:
Vu, w, a Tx (w; = (u) (w)):)),
(98)
and for the deterministic Z A case we have
Yu,a (x (w)
(99)
hence T and I produce the same counterfactuals, concluding the proof.
24

---

B.3.2 Instrumental Variable Setting
Even when an unobserved confounder is present, counterfactual identifiability can still be recovered from purely observational (C1) data, provided we can locate suitable instrumental variables. Here, an instrumental variable is any variable (or set of variables) that is: (i) statistically independent of the latent noise U: and (ii) influences the outcome exclusively through its effect on its parents PA.
Lemma B.14 (IV La-Equivalence (Nasr-Esfahang et al., 2023)), Let the instrument take values in the finite set I = {i1...., in and the endogenous parents of take values PA = {pa.,.... , Two bijective models T and T yield the same counterfactuals (i.e. are equivalent in the sense that they coincide up to a reparameterisation of the latent noise), whenever the following conditions hold:
1. (Instrumental Variable) I IL U and I 1L U.
2. For all pa E PA, the maps (:; and T' pa) are strictly monotone (either strictly increasing or strictly decreasing) and twice continuously differentiable.
3. The density (.).) continuously differentiable.
4. For every i, pa the joint density pa(i, pa. is continuously differentiable in its latent coordinate.
5. (Positivity) Vu, à and every pa PA, PU,PA(U, pa) > 0 and PŪ.PA (ú, pa) > 0.
6. (Variability) For every fixed the matrix
PD (paa (u,i))
pp(pa, (u.1))
MD(u,I)
(100)
pD u, kn)
pD(pa, u, (,
satisfies Idet (u, 1)I c, for some constant c > 0, that does not depend on u.
Proof. The proof (for d 1) is given by Nasr- Nasr-Esfahane Esfahany et al. (2023). Since and î are assumed to be strictly monotone only in the scalar setting, counterfactual identifiability is established for = 1 only. Theorem B. shows (drawing on Brenier (1991)'s Theorem) that if there exists a convex function Φ Rdd - R, where Tuu; pa, i) = Vu pa, i). then 7, are strietly monotone (rank-preserving) in the multi-dimensional case. Therefore, by Theorem B.11, we now generalise the above IV counterfactual identifiability result for d > 1, concluding the proof.
B.3.3 Backdoor Criterion Setting
Another scenario in which counterfactual identification from observational (C1) data is achievable (under stated assumptions), despite the presence of confounding, is when there exists a set of variables that satisfies the backdoor criterion (BC) w.r.t. the pair (PA, X). where PA are X's parents. That is, Z blocks all backdoor paths between PA and X, where each such path includes an arrow pointing into PA. The role of Z is account for all the spurious associations between PA and U.
Lemma B.15 (BC Lg-EEqivalence (Nasr- Fass-EEsahhany et al., 2023)). Let and T. be two distinct bijective models that induce identical joint distributions over the endogenous variables (PA, X). Then T and I are counterfactually equivalent (~<ta) if the following conditions hold:
1. (Backdoor Criterion) U TT PA Z and i 11 PA Z.
2. For every pa, the derivatives pa Jr (:: and pa det (=; both exist.
3. (Variability) For every there exist instances Z1. ...... Zd-1 such that det. M(,, 21, 2d41) 0,
(101)
where
Puiz(u (u 21)
V VuPoiz(u Zi)
M(u, 21,...,
22.11)
(102)
Pujz(1) za+1)
Vapu z(u zd41)
25

---

Remark B.16. The proof for Lemma B.15 is provided by Nasr-Esfahany et al. (2023) and already applies to the multi-dimensional setting (d 2 1), showing that bijectivity of T. F is sufficient for counterfactual equivalence when the BC is satisfied and the variability assumption holds.
Intuition. Because BC makes the change-of-variabe ating w.r.t. PA and stacking across d + 1 values of Z yields full-rank (variability) linear system whose only solution is zero: hence the link g(- : pa) (-: o T(: cannot depend on pa, so guu pa) = g()) and and T are counterfactually equivalent (~<g)) in the sense that they produce the same counterfactuals.
identity independent of the parents, differenti--
B.4 Assumptions: Summary & Plausibility
We open with Pearl and Bareinboim (2022):
"Assumptions are self-destructive in their honesty. The more explicit the assumption, the more criticism it invites. [...] Researchers therefore prefer to declare 'threats" in public and make assumptions in private."
Our work counters this practice: we intend to make our assumptions and constraints explicit so they can be better understood, challenged, and ultimately relaxed in light of new evidence or insight.
Markovian SCMs. The Markovianity assumption is standard in causality literature and is also the default assumption in many causal representation learning frameworks (Pearl, 2009; Hyvärinen et al., 2024). The assumption that there is no unobserved confounding is strong in most real- world scenarios. However, it is widely known that Markovianity alone is insufficient for counterfactual identifiability (Xia al., 2023; Bareinboim et al., 2022; Nasr- Essahany and Kiciman, 2023). Thus, further assumptions and/or restrictions on the functional class are necessary, such as monotonicity in the d = 1 case (Nast-Esfahany et al., 2023). This is the gap our work fills, particularly for d 1.
OT Regularity Assumptions. We assume that dim(X) = dim(U) 2 1 as is commonplace in flow-based generative models. This is straightforward enforce in practice by embedding lower- dimensional latents using dummy coordinates. For our strictest results, we assume quadratic cost and that the source Pt and target Px PA are absolutely continuous w.r.t. Lebesgue measure, with densities pu and P bounded above and below by positive constants on bounded convex supports.
These are standard assumptions in the literature, particularly in OT, which ensure that the measures Po and PAIPA admit well-behaved densities, guaranteeing the existence and uniqueness of an OT map (Brenier, 1991; Caffarelli, 1992). However, they may not always hold, and are often stronger than strictly necessary in practice. In practice, when working with empirical distributions, it is often advantageous to approximate them by smoothing with small continuous uniform or Gaussian noise so that the resulting distribution admits a density. Real-world data are often normalised and then modelled with smooth densities on (possibly truncated) compact supports, making these assumptions a common modelling choice. Optionally, if the densities are C1 and satisfy standard two-sided bounds on convex domain, this yields improved regularity, e.g. an interior 2. Brenier potential.
We require that there exist a velocity field v2 solving aip: + V (prz)) = 0 where Pu Pv. and Pi = Pxlpa. with Lipschitz in z and in 1.² 1] x Q) (Benamou and Brenier, 2000). Neural ODE parameterisations with spectral normalisation or weight clipping satisfy global Lipschitz bounds, and so do classical kernels with bounded derivatives (Chen et al., 201 8). Uniqueness of the associated ODE's solution is satisfied by the Lipschitz condition under -indel--nnss 's thcorem. We note that in practice the time-dependent density Pt is often fixed by choice of interpolant (c.g. linear (MeCann, 1997)). By further optimising the interpolant, it is possible to solve Eq. (47) under mild assumptions on the
density (Albergo and Handen-Eijnden. 2023).
For our most general result, which generalises scalar quantiles/ranks to d > we assume that the prior Pu be a uniform measure, e.g. (((0, 1l'). This was chosen for PIT-type reasons. However, a uniform on a d-dimennsional cube not strictly necessary for rank-preservation; other choices like the standard normal or the uniform on ball are available provided regularity conditions for OT are met (Carlier et al., 2016). In any case, weaker counterfactual equivalence relations may be obtained by relaxing the conditions under which wosnt-wise identification holds. Such relaxations are
26

---

often admissible in practice without completely sacrificing the utility of the method; this is in the spirit of many existing identifiability results (Khemakhem et al., 2021; Wasr-Esfahany al., 2023; Hyvärinen et al., 2024). Finally, the experiments we ran to validate our theory rely on Batch-OT to recover the global OT map in practice, which Pooladian et (2023); Tong al. (2024) provide convergence proof' for. However, Batch-OT is known to underperform for smaller batch sizes. Therefore, improving high-dimensional OT in terms of accuracy and efficiency remains key.
To summarise: with quadratic cost, (i) finite second moments ensure existence of an optimal coupling; (ii) the source absolutely continuous and dim(U) dim(.X), the OT is (a.e.) unique and given by Brenier's 's gradient map; (iii) two-sided density bounds on convex domains are standard, useful regularity assumptions, but are not always strictly necessary in practice; and (iv) a Lipschitz velocity field yields unique, invertible Benamou- Brennerr flow solving the continuity equation. When embedded a Markovian SCM with independent noise, this flow provides a unique, pank-preserving transport representation of the causal mechanism, generalising the notion of scalar quantiles to d > 1, and under our stated assumptions, supporting identifiability of multi-dimensional counterfactuals from observational data.
C Energy Based Models: Curl-free Flows
Proposition C.1 (Curl-free Flows.). Let (pt, Ve)eejo,1 be an admissible curve for the dynamic OT map in Lemma B.5, with U (S2; IR.) for each fixed pa RRK. Write the Helmholtz- decomposition U4 = Viha + St, where 8t is the curl component. Since (iv() = 0 by construction, the pair (pt: Vir)) satisfies dep: div(pt Vvv)) = Then, the Benamou-Brenier action A is:
2 J. (a))ll pa(u; pa) du dt 2 0,
A(p,v) - A(p, V =
(103)
hence replacing Vt by its curl-free component Vuh never increases the dynamic OT cost.
Proof. For each E [0, 1] the Helmholtz-Hodge decomposition provides unique split
J. pr(u,, pa) se(u, pa). Vn(u,, pa) du = 0, An E (2).
Ct
V t/2
St:
(104)
Orthogonality in (p:) then implies
||ve|12 = 1112=11Vaa 11
||
(105)
hence
S.// |/ (u; pa) 12 pe, pa) du
A(p A(p, Vτ) +
(106)
The integral above is non-negative and vanishes iff s; = and therefore, A(p,o) 2 A(p, Vv), with strict inequality for every con-conservative field.
Remark C.2. The unique minimiser of the dynamic OT problem is curl-free, as its time-1 map is the gradient of convex potential. Thus, the above result suggests that energy-based model (EBM) parameterisation of i (e.g. see De Sousa Ribeiro and Glocker (2025); Balcerak et al. (2025)) could prove to be a useful inductive bias for some problems, as it restricts the search space to curl-free vector fields. In practice, we operationalise this model by simply taking the sum of a neural network's 's output as the energy, then taking gradients use this model class primarily for the constructed counterfactual ellipse scenarios, as it can be quite computationally intensive for high-dimensions images, given our resources. We consider scalable curl-free flow models fertile ground for future work.
the input æ, yielding a curl-free vector field. We
27

---

D Counterfactual Soundness Axioms
When counterfactual ground truth is not available, we follow prior work (Monteiro et al., 2023; Ribeiro et al., 2023) and measure counterfactual composition, effectiveness and reversibility, which are axiomatic soundness properties of counterfactual functions that must hold true in all causal models (Halpern, 1998; Galles and Pearl, 1998; Pearl, 2009). We provide gentle introduction next for completeness; for more detailed treatment, please refer to Monteiro et al. (2023).
We point out that, while certainly useful diagnostic tools, the following metrics alone do not imply identification and should not be construed evidence of causal validity.
1. The composition axiom states that intervening on a variable to have value it would have had without the intervention should leave the system unchanged. In practice, composition is often measured using reconstruction error under null-intervention on the parents:
E Pa. || Tpa, pa.
(107)
We report this metric to enable fair comparisons with prior work. For a more stringent measure of composition that tests path independence of intervention sequences, we suggest:
pa). [T pai, Tpaaj, x)) - T(pa pa, pa.x)))
x>]h
(108)
Other composition tests could include comparing outcomes of one-by-oee one-by-one sequential inter- ventions against equivalent simultancous ones.
2. The effectiveness axiom states that intervening on variable have a specific value will cause the variable to take that value. Effectiveness is typically the most challenging axiom to obtain an unbiased measure of practice, as ideally one needs access an orucle to tell us when/if our interventions are faithful. In practice, pseudo-oracles in the form of parent predictor models are used to determine whether interventions are successful:
(a..a)>~ d(O((T (papppa"ppppppp
pa, x)), pa))
(109)
where O(-) is a pseudo-oracle function that returns the value of the parent(s) given the observation. Here dl-) is appropriate distance function dependent on the parent type. It is common practice to also apply the same metric to variables that are not descendants of the intervened-on variable, invariance there serves diagnostic measure of minimality.
3. The reversibility axiom precludes multiple solutions due to feedback loops, and follows directly from composition recursive systems such as DAGs. This can intuitively understood as measuring cycle consistency w.r.t. interventions on the parents:
E P(z,pa)~Pan
||
T(pa, pa", T(pa", pa. 2)11
(110)
Using less cluttered notation, we measure the extent to which the reversal x, matches the observation after n. > intervention cycles, where each cycle given by:
pa- Tp (x),
Tr pa
pa (x*).
(111)
Although often too costly in practice, recursive measurement of all the above metrics can yield more fine-grained insight into the faithfulness of interventions, biases of the associated counterfactual function, and incremental loss of identity w.r.t. the original observation (Monteiro et al., 2023).
28

---

E
Counterfactual Ellipse: Extra Results
E.1
Dataset Details
As outlined the main text, our ellipse dataset is built on Nasr-Esfahany et al. (2023)'s setup. We reiterate here for completeness. Let U IR² be the semi-major and -minor parameters of an ellipse, PA E (0, 2x) be an angle specifying a single point on the ellipse, and X E IR² be its cartesian coordinates. The data-generating process is defined follows:
Z = E21
C Uniform(-O0 00.5, 0.5)
PA (1.44254843 + 0.59701923 + Epa) mod (2TT),
Epa ()(
Un == exp (xp.6.699852742 0.2656131) Cuo
Gua Beta(1, 1)
U1 = (o(1 + cu1 exp (Sp...63323358 - 0.18070237))
Euy Exponential(1)
X₂ = (2 + (PA)
X1 = U/ + cos (PA)).
By construction, U TL PA Z, and Z satisfies the backdoor criterion (BC) w.r.t. PA - X. To induce Markovian setting. we simply randomise PA, yielding marginal independence U PA.
For the front-door criterion (FC) to apply, we modify the data-generating process slightly to include a mediator variable M, such that PA M X:
== Ex:
E Uniform(- 0.5)
PA (1.44254843z + 0.59701923 + Epa) mod (2T),
épa
(O
Ma == 8 (sin(PA) + Ema)
Emo 0.000.001)
Mj == 8 (cos(PA) + )m1)
fmm ,N(000001)
Uo exp(1 (16..96455274-" 0.2656131) + tuo
(1
Eug
U1 = (1 + 64
O. 0.18070237)
Guz Exponential()
Xo = (2 + Mb)
Xi (2 + M1).
where 8 1/v m3 + m projects the noised points back onto the unit circle. Note the non-linear confounding PA < 7 - U - X. This dataset is only used for the front-door criterion experiments, where we intentionally omit the unobserved confounder Z when building the model.
In all cases, we split our datasets into 70/10/20% for training, validation and testing, respectively.
E.2 Architecture and Experimental Setup
For all our model variants, we used the following simple residual multi-layer perceptron (MLP) parameterise the (conditional) vector field. We start with linear projection layer which takes as input the concatenation of æ. pa and a time index E [0, 1]. The input dimension is therefore We then have 3 of the following (residual) blocks:
dim(X) + dPA)
x ト LINEAR o SILU o LINEAR (AYERNORM(Z))
(112)
The hidden dimension of these blocks is 256. The head the MLP consists final layer norm and linear projection the output dimension dim(XX built by taking the mean of the MLP output, then taking gradients w.r.t. the input .r.
2. The energy-based flow model variants were
We used Py Torch (Paszke et al., 2019) to train all our modes for 500k steps, under identical hyper- parameter setups. These were determined with a light sweep over the learning rate and network width. We used the Adam optimizer with a learning rate of 10-1 weight decay of 10-4, B 0.9, 31 = 0.999, 10-8, and batch size 256. The best model was selected based on the counterfactual error (HAPE) the validation set during training. Final performance is reported on the test set.
29

---

E.3
Qualitative Results
G00.PA PA
0.50PAP
1,00
10.
BBM
100
-10
1.0
0.5
FLow
0.0
015
1.0
FT-F DT-EE
005
0L0
1.0
100
1.0
405
0.0
0L5
Figure 5: Visualising the curl (V x vx) of the learned vector field of different models over time (scaled to 1, 1),, for given intervention on the parents PA curl-free (irrotational) by design, the vector field is given by the gradient of scalar potential. We also see that the OT vector fields are smoother and OT-FLOW exhibits milder 'rotations' compared to FLOW. This is consistent with convergence the Brenier map, which is itself curl-free.
We can see that the EBM variants are
Obeervatica s)
True CF
Peod CF (NAIVE)
Pied CF (Ous)
KAREE 0.52
MAPE: PAPI: B.SS G8
KADe: 0.58
Parr: 46.24
6.65 $6.14
S12T
ST.UT
2
2.5
5.0
7.5
25
5.0
7.5
2.5
5.0
7.5
10
2
7.5-
10-
Rlnn: PATE:
3-
RAPSOAL
Pars Pare:
Pars: 1.36
Pare: Hare:
PAPE 56.1
5.0-
Pare: 50.09
2.
5.0
20
4-
Hare Rans 41.55
Plans: Pains 18.81
HAIE PAPS:
3
Parg Pare (8.9
Hann 0.73 RAnG 51.19
PAR Pare:
5.0-
30
in
Figure 6: Inferred counterfactual (CF) ellipses using our OT coupling flow vs. the naive OT coupling approach explained in Section 5. Not cherry-picked.
30

---

EBM (MaT:: 0,580)
FLOW (paas
OT-EBM (Pac: 0.000)
OT-FLow (Ma 0,001)
-2
EBM 0.492)
FLOW (pag:: 0.433)
OT-EBM 0.001)
OT- FLoF (MAS 0.001)
-2
-2
FBM (JAE 0.633)
FLOW (PAg 2000
OT-ERM (MAe 6.000
OT-FF.oW (MAE 0.000)
ERM (MAE 0.302)
FLoW (1200:0.003)
OT-EBM (1a0:0.000
OT-FFoow (MAA 0.000)
IBM
FLOW (P31:0.313)
OT (.000
OT-FFoW (RAE 0.000)
1313E (MANE: 0.806)
(poce 0.796)
OT-EBM C.000)
OT-FFLW (MAE 0.001)
E
G-
EBM (MAR:: 0,1)
FLOW (pag:: 0.148)
WTT-EBMAN 6.00:
OT-FLow (MAE (.C00)
-1
EBM (JoE 0.125)
FLOW (pag: 0.0g:0.1200 0.125)
OT (Pag: 0.000
OT--FLoW (RAE 0.000)
0
Obeeration (X
2, path
Reversal pach (20x)
Figure 7: Extra counterfactual reversibility results (not checry- Our OT coupling flow exhibits near-perfect counterfactual reversibility, satisfying the soundness axiom (Monteiro et al., 2023). Conversely, standard flows (and diffusion models, which equivalent (Gao et al., 2025)) exhibit comparatively poor reversibility upon repeated cycles, with u² often collapsing to a single point for any initial observation æ. In simple terms, straight paths are reversible; therefore, ought to learn straight paths induced by OT in order to satisfy axiomatic reversibility of counterfactual functions.
31

---

0.00 PA
0226 1,00
0.50 PA 1.00
0.75 PA 1.0
1.0 1.00
105
BBM
000
000
-10
PA
02..PA 100
0.50 PA 1.0
0.55 PA 1.0
D0.PA 1.00
0L5
FLow
000
-05
-10
0000 PA
FA 1.0
0.50 PA 1.00
0.50 PA L00
100, PA 1.00
OT-
0.5
0.0
0L5
1.0
PA 1.00
PA
P1
1.00
FL.F
-05
Figure 8: Visualising the curl (V vz) of the learned vector field of different models over time (scaled
to 1-1, 1)), for a given intervention on the parents PA = 1. We can see that the EBM variants are curl-free (irrotational) by design, the vector field given by the gradient of scalar potential. We also see that the OT vector fields are smoother and OT-FLOW exhibits milder 'rotations' compared FLOW. This is consistent with convergence the Brenier map, which is itself curl-free.
Obeervation (X =
True CF
Pred (EBM)
Prud (Ous)
10
0.64
0.58
0.76 1.00
0.30 PAPE
PAPE: PAPE: 1.10
0.51 PAPE:
LBB
20
2.5
s.o
5.0
G
Pazt
PArSi
PArE 1.10 HATE:
a-
Plas
PAaE: PAPE:
5.0
PAPE:
Part
jo0e::200 2.05
5.
3-
7.5
Paze O.TT
RATe 2:45
Harn: 0242
2.27 2.27
2
PAPI 5.30
PAPK 0.50
PAre: 0-44 RAPE: 2:35
Pare: 1.10
Parr: L86
PAPE 0.32 PAPs: 2.58
5.0-
10
2.5
5.0
Figure 9: Inferred counterfactual (CF) ellipses using our (Markovian) OT coupling flow vs. an energy-based curr-frree) parameterisation of flow- model. Not cherry-picked.
32

---

0.00 PA
0.25 PA
0.50
0.00
PA 0.00
EBM
1.0
FLow
1001
OT- ET-EBM
FLoW
Figure 10: Visualising the curl (V x U2) of the learned vector field of different models over time (scaled to [-1, 1)). for given intervention on the parents PA = -1. We can see that the EBM variants are curl-free (irrotational) by design, as the vector field is given by the gradient potential. We also see that the OT vector fields are smoother and OT-FLOw exhibits milder 'rotations' compared to FLOW. This consistent with convergence to the Brenier map, which itself curl-free.
scalar
Obeervation (X =
True
Prod CF (FL.ww)
Pred CF (Gus)
Hare: Hare:
Alare: MAme:
HAPE: PAPE:
PaPS: PAPS
PAric Pare:
PArE:0.51 PaPe:2.00
2
dn
5.0
7.5
7.5
Pare
Para 1.51 Pars
Pare:: pares
Pare: 2.P::2.25
PATE pares 3.87
Pare
5.0
Rope:
5.0
2.5
3-
7.5
Pasep O.77 (4400: 3.07
Hapy:: 1t 2.00::205 2.85
PAPR (1522 MAPRIA.7
Rape: 0.44 (AAe 2.50
pape: Pare:111 1.10 3.P
Pars: 11R2 PAPE: 8.30
10
10
Figure 11: Inferred counterfactual (CF) ellipses using our (Markovian) OT coupling flow vs. continuous-time flow (Lipman al., 2023). We note that standard continuous-tine flow equivalent to diffusional model (Gao et al., 2025); thus, this baseline is also representative of Sanchez and Tsaftaris (2021)'s 's diffusion-bart based counterfactual inference model. Not aherry-pickeed
33

---

E,4
Ablation Studies: Priors & Assumption Violations
Table 3: Ablation of complex data-generating processes (DGP) CF error μAPE (%)
Counterfactual Ellipse
DGP Pt
EBM
FLOw
OT-EBM
OT-FLOW
Original
2.32
2.30+.02
0.93#.02
0.76+.01
Bimodal
5.112.03
5.01
1.73#.09
1.49
Multimodal
7.94
7.69
2.29#.02
2.03#.03
More Complex Pv. Intuitively, more complex data-generating process is expected be harder to model. To test this hypothesis, we created new datasets with multimodal Pv- We use an indicator 8; Bcrnoulli(055)
induce non-linear, bimodal shift into μ4; E (-2, 2} into U's mechanism. We also use a categorical shift μu; E (-4. 2, 2, 4}, with respective shift probabilities (0.3, 0.2, 0.2, 0.3). generate another multimodal setting. As reported in Table 3, we confirm that more complex Pt can affect generation, but importantly, the relative rank of our models remains consistent, Naturally, we expect larger models to reduce performance gaps for complicated Pt's.
Table 4: Ablation of monotonicity violation.
Counterfactual Ellipse
MODEL
Monotone
COMPOSITION | REVERSIBILITY 1
НАPE (%)
Nasr- Basr-Esfahany al. (2023)
N
607.NIA
EBM
N
27.86+ 369
39.15+.349
2.319
FLOW
27.99+.08 0096
39.39+.229
2.295=
OT-EBM
Y
0.659
1.056=
0.925=.019
OT-FLOW
0.618:1.066
1.016:1.000
0.763
Monotonicity Violation. We ran new set of experiments to evaluate different methods under a monotonicity violation. In addition to measuring counterfactual error, inspect composition and reversibility counterfactual soundness axioms, in order to reveal possible performance disparities. For evaluation, reversibility measurements. As reported in Table 4, models which are not "vector-monoton" substantially inferior counterfactual soundness. The results confirm the intuition that straighter paths lead to improved reversal capabilities, That said, whether this generalises to much higher dimensions entirely dependent the quality of the dynamic OT approximation being used.
50/250 steps for ODE/SDE solving, and 20 cycles for composition and exhibit
Table 5: Ablation of bijectivity violation.
Counterfactual Ellipse
MODEL
Bijective
COMPOSITION |
REVERSIBILITY &
РAPE (%)
SDE Abduction
N
2.9504.013
4.173+.019
3.0644.011
SDE Prediction
N
3.006#.005
4.276#.019
3.126#.000
CFG 1.05)
N
0.032+.001
0.487#.010
2.174+.0m
CFG (w = 1.15)
N
0.032+.001
1,499 0013
6.081+.000
CFG (w 1.25)
N
0.032
2.598
9.809
OT-FL.W
Y
0.032.1.04
0.05611.09
0.763.1.014
Bijectivity Violation. In this experiment, we use three strategies to stress-test bijectivity violations: (i) stochastic abduction by solving an SDE instead of an ODE; (ii) stochastic prediction by solving an SDE; (iii) classifier-fres free guidance (CFG) (Ho and Salimans, 2022). One cycle was used for measuring composition and reversibility. As reported in Table 5, bijectivity violations affect reversibility and composition substantially; this is likely due to the ground-truth mechanism being bijective.
34

---

Table 6: Ablation of Markovianity violation. The NAIVE OT coupling baseline corresponds to Batch-OT flow matching (Pooladian et al., 2023; Tong et al., 2024).
Counterfactual Ellipse
OT Coupling
Markovian
COMPOSITION + REVERSIBILITY + μAPE. (%) +
NAIVE
N
0.063
0.346#.012
4..01
Ours
Y
0.0322.004
0.056#.009
0.763+.015
Markovianity Violation. In this experiment, we demonstrate the importance of using the correct model specification that reflects causal assumptions. As explained in the main text, the naive Batch-OT coupling violates the Markovianity assumption, so it should not be used for counterfactual inference if the assumption were hold. As reported in Table 6, observe that composition and reversibility are reasonable for the NAIVE OT coupling version, as the map remains bijective. How- ever, the counterfactual error is substantially higher for the NAIVE version compared to our proposed Markovian OT coupling, since the Markovianity violation leads to inconsistent counterfactuals.
F
MIMIC Chest X-ray: Extra Results
F.1
Dataset Details
We reproduce the dataset setup used by our baselines, i.e. Ribeiro et al. (2023) and Xia et al. (2024). The chest X-ray images were resized to 192 x 192 resolution. From the associated metadata, we focused on four key attributes: SEX (S), RACE (/t), AGE (A) and DISEASE (D). The assumed causal graph follows prior work (Ribeiro et al., 2023) in that A + D and PA = (S, R. A, D} are the causal parents of the chest X-ray image X. We also focus our scope to only pleural effusion for the disease, to keep comparisons fair. As result, the final dataset includes only individuals who were either diagnoses with pleural effusion (diseased) reported having findings (i.e. healthy). We then divided the dataset into 62,336 subjects for training, 9,968 for validation and 30,535 for testing, again following the exact same protocol as both Ribeiro al. (2023) and Xia et al. (2024).
F.2 Architecture and Experimental Setup
The causal mechanisms for all variables PA R, R, A, D3 and X were learned from observed data. For PA, standard discrete-time time flow-based modelling was used following Pawlowski et al. (2020); Ribeiro et al. (2023). For X, we use scaled-up version of our most successful model analysed the counterfactual ellipse experiments, namely OT-FLOW, which includes our specialised family of batch TO couplings to satisfy the Markovianity requirement (ef. Section 5). With that said, we recognise that using large batch sizes is best for accurate batch OT approximations in high dimensions (Klein et al., 2025). Given current resource constraints, we encourage future work to explore more accurate and efficient OT approximations at scale, possibly along the lines of Mousavi-Hosseee et al. (2025).
All models are trained in continuous time using our Markovian Batch-OT coupling. Since the AGE (A) attribute is continuous, we use an age binning strategy to ensure we sample from the conditional distribution PC P as explained in Section 5. The effect is that each sample in each batch has the same parents, except for age, which may deviate slightly depending on the bin size. Table 12 reports an ablation study on the age bin size; we find negligible performance differences across bin sizes.
To parameterise the (conditional) vector field for all our models and ablation study variants, we used the same streamlined version of the UNet proposed by Dhariwal and Nichol (2021). To condition the model on the parent variables, we simply learn a separate embedding vector for each parent and sum them all together. This joint parent embedding is then combined with the time embedding to condition each block in the UNet, following Dhariwal and Nichol (2021). The hyperparameters are given in Table 7. The model has just 22M trainable parameters. For training, we used the AdamW optimiser with learning rate of 10-4, weight decay of 10-4, 33 = 0.9, 31 = 0.999, 10-*, and a batch size of 64. We use a linear learning rate warmup of 2000 steps. No extensive hyperparameter sweep was necessary to obtain sufficiently good performance. We trained our models for a maximum of 300 epochs. Since our evaluation is focused on the counterfactual soundness axioms, which are
35

---

quite costly to evaluate, we instead performed model selection using the validation set loss achieved by an exponential moving average (EMA) of the model parameters (EMA rate of 0.9999). Further improvements in performance are expected from using standard sample quality metrics for model selection. We conducted three runs with different random seeds. All our models were trained on L.40 GPUs, with the full model fitting entirely on a single GPU.
F.3 Additional Comparisons and Ablation Studies
Tables 8 and report additional comparative results against the baselines (Ribeiro al., 2023; Xia et al., 2024). We again find that our method is superior for all three measured counterfactual soundness axioms, often by large margins, without requiring any costly counterfactual fine-tuning or classifier(-free guidance strategies. We also conduct ablation studies on our proposed modifications, and the results are reported in Tables 10, 11 and 12. In summary, we observe significant improvements from using our OT coupling compared to the NAIVE non-Markovian approach. As outlined in the main text, we find different performance trade-off's between OT-FLOW and FLOW for different interventions, suggesting either Markovianity violation subpar OT approximation (due to small batch size of 64). We expect further improvements from increasing model capacity and batch size significantly (Klein et al., 2025), to better approximate the global OT map (Pooladian al., 2023).
Table 7: UNet hyperparameters.
MIMIC Chest X-ray
image shape
(1, 192, 192)
model channels
32
channel mult
[1, 4, 6, 8]
residual blocks
2
attention resolutions
[-1]
heads
1
head channels
64
dropout
0.1
Table 8: Measuring counterfactual effectiveness on MIMIC Chest X-ray. |Auc denotes the absolute difference in ROCAUC of inferred counterfactuals relative to the observed data baseline. For each variable, our results (blue shade) appear on the right, and baseline results are on the left. Unlike the baselines, which rely on costly fine-tuning stage with classifiers/regresoor for each variable, our approach does not require any form classifier or classifier- free guidance to perform well.
MIMIC Chest X-ray (192x 192)
SEX (S)
RACE (R)
AGE (A)
DISEASE (D)
BASELINE
AUC (%)
AUC (%) f
MAE (yr) 4
AUC (%) 1
Observed data
99.63
95.34
6.197
94.41
(w/o CF fine-tuning Ribciro et al. (2023))
INTERVENTION
AAuu| (%) $
Axve (%) 4
Amar (yr)
Aave| (%) |
SS 8) 0.130
7,430
0.173
20.84
0.583
0.486
0.333
1.310
0.023
do R = r)
0.180
36.94
0.050
0.229
0.394
7.810
0.310
do(A A = a)
0.130
0.
20.44
1.197
3.872
0.836=.08
6.110
0.347
do(D D d)
0.030
0.067 4.00
20.14
0.627+.18
0.560
0,435 s..5
22.01
2.280
do(rand)
1.330
0.150+.02
24.14
0.730
1.166
0.510 a..5
7.010
0.640+..1
(Soft CF finc-tuning Xia al. (2024))
INTERVENTION
AAuc (%) 4
|Axeel (%) 4
AMAL (yr)
A (%) 4
SS = 8)
0.070
0.173 oz
2,040
0.583
0.333
0.090
0.023
do R = =
0.130
0.180
3.360
0.050 07
0.394
0.110
0.310
do(D- = d)
0.030
0.066.1.0
1.540
0.627.1.18
0.435 1.05
3,690
2.280
36

---

Table 9: Measuring counterfactual reversibility and composition properties (Monteiro et al., 2023). The baseline results are reproduced exactly according to Ribeiro et (2023). We use an Fuler ODE solver or dropi5 (le-5 tol). A single reverse cycle was used for computing reversibility. Three randomly seeded subsets 1000 test samples were used, and the results were averaged. Similar initial performance was observed for FLOW, thus improving OT further is expected to boost performance.
MIMIC Chest X-ray (192x 192)
METHOD
NFE
COMPOSITION MAE (px) 4
REVERSIBILITY MAE (px) +
HVAE (Ribeiro et al., 2023)
N/A
3.09543 1 0.0536
2.89816 0.4153
2,70379 + 0.1151
20.6228 0.6492
OT-FLoW (Ours)
250
0.55289 1 0.0280
1.22693 0.1630
dopri5
0.18352 + 0.0276
0.49485 + 0.1549
Table 10: Ablation study of our conditional OT coupling counterfactual effectiveness. The baseline is the standard OT coupling (NAIVE) described Section 5, which violates the Markovian requirement. The results show approach boosts counterfactual effectiveness, especially for disease and race interventions. The most notable improvements are highlighted using shaded blue areas.
MIMIC Chest X-ray (192x 192)
SEX (S)
RACE (R)
AGE (4)
DISEASE (D)
INTERVENTION
OT Coupling
AUC (%) t
AUC (%)
MAE (yr) 4
AUC (%)
N/A (Observed data)
N/A
99.63
95.34
6.197
94.41
NAIVE
99.42
93.85
6.602
94.47
do(S = s)
Ours
99.45
94.761
6.529
94.39
NAIVE
99.56
88,03
6.482
94.24
do(R ==r) =r)
Ours
99.45
95.29
6.590s.10
94.10
NAIVE
99.53
93.95
7.723
93.88
do(A = a)
Ours
99.442.01
94. 1.4 16
7.0322
94.762209
NAIVE
99.58
93,66
6.630
88.35
do(D d)
Ours
99.53=
94.71
6.6322.05
92.13
NAIVE Ours
99.55
92.75
6.844
92.62
do(rand)
99,48-
94.61+:23
6.706+
93.77+t117
37

---

Table 11: Ablation study of Batch-OT on counterfactual effectiveness. We observe different per- formance trade-offs for different models; e.g., the OT-FLOW performs better for race interventions while performing worse than FLOW for disease, suggesting possible non-Markovian interactions. We note that batch size of 64 was used for training due current GPU memory constraints. Larger performance differences are expected with larger batch sizes (Klein et al., 2025), or by adopting more accurate and scalable OT approximations Mlousavi-Hosseeee et al., 2025).
MIMIC Chest X-ray (192x 192)
SEX (5)
RACE (R)
AGE (A)
DISEASE (D) AUC (%) 1
INTERVENTION
OT
AUC (%) f
AUC (%) f
MAE (yr) +
N/A (Observed data)
N/A
99.63
95.34
6.197
94.41
99.58
94.64
6.495
94.72
do(S s)
Y
99.45
94.762.15
6.529 06
94,39
99.56
93.09
6,460
94.59
do(R=
Y
.9
95.29
6.590
94.10
N
99.54
94,40
6.855
95.23
do(A a)
99.44
94.14..16
7.032:.07
94.76.1.9
N
99.56
94.26
6.589
92.95
do(D d)
99.53
94.71
6.632
92.13
N
99.56
94.17
6.596
94.48
do(rand)
99.48
94.61
6.706
93.77
Table 12: Ablation study of AGE (A) bin widths on counterfactual effectiveness. We observe minimal differences performance across different bin sizes. This is encouraging as it demonstrates our OT coupling strategy can still work well for continuous parent variables.
MIMIC Chest X-ray (192 192)
SEX (S)
RACE (R)
AGE (A)
DISEASE (D)
INTERVENTION
AUC (%) ?
AUC (%)
MAE (yr)
AUC (%)
N/A (Observed data)
N/A
99.63
95.34
6.197
94.41
99.48
94.49
6.580+
94.37+
do(S =
99.45
94.76
6.529
94.39
5
99.41+.02
94.85
6.600+
94,44
05
94.49=1.1
6.6054 o5
94.07:
do(R =r)
99. 05
95.29
6.5904 10
94.10=:19
92
94.97
6.625 od
94.03:1.0
2
a9
94.06
7.1374.00 .6
94.674
do(A
99.44
94.14
7.032:1.0 07
94,762.09
99.43
93.89
7.159 05
94.62:1.02
99.56++.zz
94.44
6.675
91.72+
do(D = d)
14.71=.19
6.6324 os
92.133.37
99.52+.m
94.70=223
6.676+ 03
91.664
99.49
94.35.1..0
6.7554.00
93.774:24
do(rand)
3
94.61
6.7064 04
93.77+.17
99.46 5
94.53
6.760 o
93.664
38

---

F4
Qualitative Results
Dlratlay
Ohoess Obeesration
MeSS - Tanak?
dil badital
de.5 - (ando))
wit8 - caie)
dee mins
dee mn diomall
AolA m Hkx
Coucseefacke tanl
anla aelo
Reversa
[Twckle rode wike)
ckite
(wallly
Cloerraiio
Ar8 Fanaii
AuS mak
Sh Initle)
D - aach
de' - fenelo
docl de)
Connater Sacteaactto
baitted
Reverall
(Joatly 2he male and Caialie, diy fmads Maak di fhey lneal wit Lnoosra oi day offeifelionll fhu Dmale all shooee 2a, ou lo skitff lovisl Sky fendl wline Bookle 70s Smake alass
Clees vation
do dewsad
Trenkel
Conssesfactun
4s10 balike
de - vtatol
de- - diol
de decuol
4:3 sabo
#i15 Fmaal
Mr Ianak
Heveral
Figure 12: Qualitative counterfactual inference and reversibility results (500 ODE Euler solver steps) using our Markovian OT-FLOW model. Despite the small model size, we observe faithful, reasonably eieetity-preserving interventions, well impressive counterfactual reversibility. Importantly, no costly counterfactual fine-tuning (Ribeiro et al., 2023), or classifier(-free) strategies were required.
39

---

[inoonn oity male vrand Donilte des makk weeid Donalles die mede wel (ooil stt lime liek (eteooeni, lit mak, ever)) Chvaltts dhr, the welend (boaityy dor mark win fohn and the keek elend
Obeetvation
Countesfacuu
O in danesi
Ald Inablyy
deS teak
de teable)
de reaai
Fonerall
eies
Jmakl baccal socnece Cmeliio
tay
Clecrvatio
do.2
dee i loaaki
dr ducaoh
46 D docsed)
Counterfacuu
6d8 Sessa)
makk
beakton
de beitter
Hecemal
Obaervaiio
0eS mavi
du beatten
Ai diwwwdii
Coennteffcctn
Heveessl
Figure 13: Qualitative counterfactual inference and reversibility results (500 ODE Euler solver steps) using our Markovian OT-FLOW model. Despite the small model size, we observe faithful, reasonably ddentity-preservie interventions, well impressive counterfactual reversibility. Importantly, no costly counterfactual fine-tuning (Ribeiro et al., 2023), or classifier strategies were required.
40

---

Observatin
detrak
Andaoweli
Aktrdd
Countrefac Coutteeeacttuu
Direct Elleet
kniityy calk wllls
lnbiiie dily wakk velio
Tealfie foods
Obeernaioon
asicuso)
(eicbounc)
tichnatte
dejoraci
Direst Counterfactos
Lasped Seado sclit
bastto calk watt
Leahta Dhicnalo stlii
beatt
Meck
luchte Donaiaa walte dissed de eaceett
Oboer votion
Countrfictuul
Direct Eiscet
Figure 14: Extra qualitative counterfactual inference results with direct effect maps (500 ODE Euler solver steps). We observe fairly good identity- preservtonnn overall, but counterfactual effectiveness appears to be somewhat prioritised by the model. We anticipate that identity preservation can be further improved using inference time guidance strategies without having retrain the model.
41

---

desand cieaal lukite dasused
Balitt Ca
The
desand rude
Obseration
detheithe)
dencavi
deraak
Counterrictuu
Direct Eliees
Twckile Thy, in wtitnn dieisel lim fmman aldss leally, dyy los dl Mall dicciel Thy linah whine dioond Toa Dmaale winte beilil,, doy main utlis Featly tht fieak whits dioonvd Ro, recls vila
Oloeetration
delliowell
dalabile)
dvllioelli
deckeade
L'xees Eibeet. Countrsfactaal
Observtiin
Aetrate)
wann
Ardwowdd
Countieffee facctttttu
Direct Aileet
Qualitative counterfactual inference results with only 50 ODE Fuler solver steps. We observe decent results considering the small number of ODE solving steps used (classical diffusion models would use, e.g. 1000). However, we find that identity preservation is more challenging in this regime, and spurious associations (e.g. background, artefacts (Pérez-Garcfa al., 2025)) tend to be more prevalent. We expect inference-time classifier(-free) guidance to further improve the results.
Figure 15:
42

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: We claim counterfactual identifiability for multi-dimensiona. variables in Markovian SCMs using our prescribed model, and validate our theoretical result experi- mentally using a carefully constructed scenario.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results be expected to generalize other settings.
is fine include aspirational goals as motivation long is clear that these goals are not attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: We urge practitioners to consider the ethical implications of their causal assumptions when applying this framework in safety-critical settings.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "I "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and prool's
Question: For each theoretical result, does the paper provide the full of assumptions and complete (and correct) proof?
43

---

Answer: [Yes]
Justification: We provide detailed appendix with both motivating examples and proofs for all our theoretical claims.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification: We provide all the details required reproduce our results, both: (i) methodology-wist in the main manuscript (Sections 4, 5 and 6); and (ii) all implementation details terms of neural architecture and hyperparameters are provided in the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how replicate the results, access a hosted model (e.g.. in the of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily
model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers to have some path to reproducing verifying the results.
44

---

5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, described in supplemental material?
Answer: [No]
Justification: We only publicly available datasets (Ellipse and MIMIC- CXR). We
working on an open-source version the code and guarantee that it will be made publicly available.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips cc/
for more details.
While we encourage the release of code and data, we understand that this might be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https:
/ /nips.cc/ppll wtt.c//pcpcgg..eeeeeeeeooeeepeee (ooeeSubiissii soimiiiissiiioiiooooccoc.mmmmmmc for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: We provide the primary details in the experiments section, and the remaining details can be found in the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment statistical significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification: We ran experiments with multiple random seeds and report the mean and standard deviation of the results.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emoor bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
45

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: GPU details are provided in the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics https: /ttee////eeueeeeeeee /nerrpss....p// cc/
Answer: [Yes]
Justification:
We have taken steps to ensure that this is the case.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- cration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: We warn users to carefully consider the impact of their technical and causal assumptions when drawing causal conclusions from our model in real-world settings.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
46

---

If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), faimess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g.. gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [No]
Justification: We only use publicly available data, and our models are small-scale by modern standards. We don't 't believe they present additional risks relative to the systems already available or warrant specialised safeguards.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g.. code, data, models), used in the paper. properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: We cite all relevant sources and only make use of publicly available academic material and datasets.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
47

---

For scraped data from a particular source (e.g., website). the copyright and terms of service of that souree should be provided.
assets are released, the license, copyright information, and terms use in the package should be provided. For popular datasets, paperswi thcode. com comddatases has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
No new assets are provided at this time.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/mmodll as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
We did not conduct such studies in our work.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals equivalent for research with human subjects
Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Question:
Answer: [NA]
Justification:
We did conduct such studies in our work.
Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with human subjects.
48

---

Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if' applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: LLMs were occasionally used for help with formatting and writing plotting code.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non-standard components.
Please refer to our LLM policy (https: / bttppa//neurre htt.iipp.cccoooffeeeecocccm ec//Oonferencee2225//M.. for what should or should not be described.
49


Reviews:
[{"review_id": "ehVIAovGHc", "reviewer": "Reviewer_hTdz", "summary": "Counterfactual identifiability ($\\mathcal{L}_3$) from observational data is both important and challenging. According to Pearl's Causal Ladder, it is generally impossible without (often strong) assumptions. Previous methods for counterfactual identifiability often require model assumptions such as linearity, bijectivity, Lipschitz continuity, or the Markovian property. They also rely on data assumptions such as faithfulness, atomic interventions, interventional data, or paired data (e.g., from randomized controlled trials or counterfactual datasets) obtained under different interventional distributions.\n\nThis paper, primarily extending the work of Nasr-Esfahany et al. (2023), builds theoretical guarantees based on bijectivity and monotonicity assumptions to extend identifiability from one-dimensional endogenous variables $X_i$ to $d$-dimensional endogenous variables $X_i$. The goal is to infer counterfactual queries of $X_i$ given its parents $\\textbf{pa}_i$.\n\nThe theortical contributions of the paper follow this roadmap:\n\n1. Proposition 4.3 shows that if the mechanism is monotone in $u$, then the counterfactual transport map $T$ is monotone in $x$, linking the data-generating process assumption to a learnable property of the function.\n\n2. Proposition 4.5 shows that, under several assumptions, monotonicity leads to a unique mapping $T$ determined by $P_U$ and $P_{X|\\textbf{PA}}^\\mathfrak{C}$.\n\n3. Proposition 4.8 attempts to prove that two transport maps are counterfactually equivalent if there exists an invertible map $g$.\n\n4. Proposition 4.10 shows the existence of a monotone, bijective map $g^*$ that minimizes the quadratic transport loss.\n\n5. Proposition 4.11, building on the previous propositions, establishes counterfactual identifiability for Markovian Structural Causal Models (SCMs).\n\nFurthermore, the paper proposes a practical matching method using dynamic optimal transport (OT) flow and demonstrates its effectiveness on both synthetic datasets and a real-world X-ray dataset.", "questions": "1. In Proposition 4.8, two transport maps $T^{(1)}$ and $T^{(2)}$ are equivalent if there exists an invertible transition map $g$, and in Proposition 4.10, two prior distributions can be mapped by a unique, bijective, and monotone transport map $g^\\*$ that minimizes the quadratic transport cost. But in Proposition 4.11, we need to identify/learn the $g$ between the two priors: the uniformly distributed prior and the unknown true prior. We cannot directly learn the exogenous prior from data; thus, minimizing the quadratic transport cost might lead to a $g$ that is not the true $g^\\*$ we want.\n\n2. What is the assumption on the interventions: soft intervention or hard intervention? It seems that the paper requires the intervention to preserve the monotonicity of the interventional SCM. Is that too restrictive?\n\n3. In the synthetic experiment, what is $Z$? (line 287)\n\n4. Intuitively, in section 6.1, given a point observation X=x, how possible could we learn the whole ellipse? Could you elaborate more on the experiment and connect it to the theoretical assumptions and guarantees?", "limitations": "this paper discusses the limitation briefly in the conclusion section.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "**Strengths**\n\n1. The extension of counterfactual identifiability in Bijective Causal Models from one-dimensional to $d$-dimensional endogenous variables is a significant and meaningful contribution.\n\n2. The paper is clearly written, with minimal typographical errors.\n\n3. The paper is well-organized. In particular, Section 4 presents a series of theoretical results in a clear and structured manner.\n\n**Weaknesses**\n1. The strong assumptions and setting (e.g., known DAG, bijective and monotonic SCM) might limit the practical usefulness of the method. Could the authors elaborate on potential application scenarios where these assumptions are reasonable?\n\n2. My main concern for the synthetic dataset is about the selection of $P(U).$ If a more complicated $P(U)$ is used in the data generating process, would the generation result be affected?\n3. Due to the restrictive assumptions made in this paper, additional experiments might be beneficial. For example, how does the proposed method perform when certain assumptions are violated? Specifically, without the monotonicity assumption or without bijectiveness, to what extent is performance affected?\n\n**Minor Comments**\n\n1. In Definition 4.1 (Markovian SCM.), the period inside the parentheses should be removed.\n\n2. On page 3, line 106, what does $P_X^{\\mathfrak{C}_{x'}|X=x}$ mean? Particular, how is an SCM conditioned on $X = x$ interpreted? Clarification would be helpful.", "was_revised": true, "final_justification": "Counterfactual identifiability from observational data is both important and challenging. By extending the analysis from 1-dimensional to a $d$-dimensional endogenous variable $X_i$, the authors take this challenge one step further. I recommend accept for this paper due to its theoretical contribution. Notice that I didn't check the proofs in the appendix.", "rebuttal": "We thank the reviewer for their positive feedback. We are pleased that our counterfactual identifiability result was deemed **\"a significant and meaningful contribution\"**. We also appreciate the recognition of our paper's clarity, minimal typographical issues, and well-organised structure, particularly the presentation of theoretical results. \n\n> ### **Q1: The strong assumptions (e.g., known DAG, bijective and monotonic SCM) might limit the practical usefulness of the method. Could the authors elaborate on potential application scenarios where these assumptions are reasonable?**\n\nIn **Section 6.2**, we demonstrate the practical utility of the method on a high-dimensional medical imaging problem, outperforming existing works by a large margin.\n\n**Known DAG:** Many causal relationships are known in the medical domain, thanks to a long history of results from randomised controlled trials and meta-analyses. This positions healthcare as the ideal application area for causal modelling, as one can tap into existing medical knowledge to help build causal graphs.\n\n**Bijectivity:** Bijective causal mechanisms (Pawlowski et al. 2020, Nasr-Esfahany et al. 2023) subsume many model classes with known identifiability results, e.g., (non)-linear additive noise models (Shimizu et al. 2006; Hoyer et al. 2008; Peters et al. 2014), post-nonlinear models (K. Zhang and Hyvärinen 2009), and location-scale models (Immer et al. 2023). As such, they represent an attractive model class with practical instantiations using modern flows.\n\n**Monotone Map:** Learning high-dimensional monotone maps is challenging. However, this is now provably possible in practice thanks to recent advances (Pooladian et al. 2023; Tong et al. 2024). When the model is misspecified and/or the assumptions are violated, guarantees may not hold. This is not unique to our approach, but holds for all identification results.\n\nMany recent works already rely on such assumptions (Pawlowski et al. 2020; Sanchez and Tsaftaris 2021; Monteiro et al. 2023; Ribeiro et al. 2023), but their causal claims remain speculative without identification, representing an operational risk. Thus, our work represents an important advance for the responsible use of causal AI.\n\nWe also provide theoretical results for non-Markovian settings (**Section 4.3**), when 3 common criteria apply: Backdoor Criterion, Instrument Variable and the Frontdoor Criterion, rendering our approach widely applicable.\n\n> ### **Q2: My main concern for the synthetic dataset is about the selection of P(U). If a more complicated P(U) is used in the data-generating process, would the generation result be affected?**\n\nIntuitively, a more complex data-generating process (DGP) is expected to be harder to model. \n\nTo test this hypothesis, we created new datasets with multimodal $P_U$. We use an indicator $s_i \\sim \\text{Bern}(0.5)$ to induce a non-linear, bimodal shift $\\mu_i \\in \\\\{-2, 2\\\\}$ into $U$'s mechanism. We also use a categorical shift, with $\\mu_i \\in \\\\{-2, 0, 2\\\\}$, to generate another multimodal setting. More details on all new experiments will be added to the paper.\n\nWe confirm that more complex $P_U$'s can affect generation, but importantly, the relative rank of our models remains consistent. Naturally, we expect larger models to reduce performance gaps for complicated $P_U$'s.\n\n|DGP $P_U$|EBM|Flow|OT-EBM|OT-Flow|\n|-|-|-|-|-|\n|Original|2.32±.01|2.30±.02|0.93±.02|0.76±.01|\n|Bimodal|4.95±.01|4.95±.08|1.72±.05|1.40±.05|\n|Multimodal|5.11±.03|5.01±.03|1.73±.09|1.49±.06|\n\n> ### **Q3: How does the proposed method perform when certain assumptions are violated? Specifically, without the monotonicity assumption or without bijectiveness, to what extent is performance affected?**\n\nTo thoroughly address this, we ran a new set of experiments.\n\nIn addition to measuring counterfactual error, we inspect the **Composition ($\\downarrow$)** and **Reversibility ($\\downarrow$)** counterfactual soundness axioms, to reveal performance disparities under assumption violations. \n\nWe use the Mean Average Percentage Error (MAPE) metric, and 50/250 steps for ODE/SDE solving.\n\n**Experiment 1: Monotonicity violation**\n\nWe use 20 cycles for composition and reversibility measurements; models which are not vector-monotone exhibit significantly inferior counterfactual soundness.\n\n|Model|Monotone|Composition|Reversibility|CF Error|\n|-|-|-|-|-|\n|Nasr-Esfahany et al. (2023)|N|-|-|607|\n|EBM|N|27.86±.369|39.15±.349|2.319±.006|\n|Flow|N|27.99±.096|39.39±.229|2.295±.020|\n|OT-EBM|Y|0.659±.020|1.056±.072 |0.925±.019|\n|OT-Flow|Y|**0.618**±.046|**1.016**±.090 |**0.763**±.015|\n\n**Experiment 2: Bijectivity violation**\n\nWe use 3 strategies to stress-test bijectivity violations: **(i)** stochastic *abduction* by solving an SDE instead of an ODE; **(ii)** stochastic *prediction* by solving an SDE; **(iii)** classifier-free guidance (CFG). One cycle was used for composition and reversibility. \n\nWe observe that bijectivity violations affect reversibility and composition significantly, as it becomes harder to 'undo' interventions.\n\n|Model|Bijective|Composition|Reversibility|CF Error|\n|-|-|-|-|-|\n|SDE Abduction|N|2.950±.013|4.173±.019|3.064±.011|\n|SDE Prediction|N|3.006±.015|4.276±.019|3.126±.007|\n|CFG (w=1.05)|N|0.032±.004|0.487±.010|2.174±.070|\n|CFG (w=1.15)|N|0.032±.004|1.499±.013|6.081±.069|\n|CFG (w=1.25)|N|0.032±.004|2.598±.024|9.809±.078|\n|OT-Flow|Y|0.032±.004|**0.056**±.009|**0.763**±.015|\n\n**Experiment 3: Markovianity violation**\n\nWe observe that composition and reversibility are reasonable for the Naive OT-Flow version, as the map remains bijective and monotone. However, the counterfactual error is significantly higher for the Naive version, as the violation leads to inconsistent counterfactuals.\n\n|OT-Flow|Markovian|Composition|Reversibility|CF Error|\n|-|-|-|-|-|\n|Naive|N|0.063±.005|0.346±.012|47.01±.016|\n|Ours|Y|**0.032**±.004|**0.056**±.009 |**0.763**±.015|\n\n> ### **Q4: In Proposition 4.11, we need to identify/learn the g between the uniform prior and the true prior. We cannot directly learn the exogenous prior from data; thus, minimizing the transport cost might lead to a g that is not the true g we want.**\n\nThank you for the thoughtful question. To clarify, learning the prior transition map $g$ explicitly is not strictly necessary. \n\nLet the true generator be $x = f(u;\\operatorname{pa})$, $u \\sim P_U$. We train a single OT-Flow $T(u;\\operatorname{pa})$, $u \\sim U(0,1)^d$, to match the output law.\n\nAlthough there exists a unique OT solution for $T$, we make no claims that it factorises exactly into $f \\circ g$, for any choice of $P_U$. There could be a map with a lower transport cost that only implicitly encodes $g$.\n\nTo illustrate this simply, consider a rubber-cube representing the uniform prior $U(0,1)^d$.\n\n1. A warp $g$ stretches the cube so that the embedded points follow the true prior $P_U$;\n\n2. A map $f$ then moulds the cube into data space.\n\nA single map $T$ could enable steps 1 & 2 in unison if this reduces transport cost. In analogous terms, a path A-to-C may be closer than a path A-to-B-to-C. We will make this point clearer in the final version.\n\n> ### **Q5: What is the assumption on the interventions: soft intervention or hard intervention? It seems that the paper requires the intervention to preserve the monotonicity of the interventional SCM. Is that too restrictive?**\n\nWe focus on hard interventions as they enable many useful real-world applications. This is exemplified by our experiments in Section 6.2, and ongoing work in: (i) evaluating the impact of interventions (Kusner et al. 36 2017; Tsirtsis and Rodriguez 2023), (ii) understanding cause-effect relationships (Karimi et al. 2020; Budhathoki et al. 2022), (iii) generating targeted synthetic data (Pitis et al. 2022; Roschewitz et al. 2024), and (iv) quantifying disease severity (Mehta et al. 2025).\n\nWhile further theoretical work is needed to accommodate *arbitrary* soft interventions, as they can alter the underlying mechanism, the current method supports a suitably defined subset of soft interventions that preserve vector-monotonicity.\n\n> ### **Q6: In the synthetic experiment, what is Z? (line 287)**\n\nWe thank the reviewer, as their question helped us find a typo. The random variable $Z$ is part of the data-generating process, a confounder satisfying the Backdoor Criterion (BC) w.r.t. $\\operatorname{PA} \\to X$. We include it to reproduce Nasr-Esfahany et al. (2023)’s BC result. In the text, we inadvertently state that we randomised $Z$ to induce a Markovian setting, but what we did was randomise $\\operatorname{PA}$, ensuring $U  \\perp\\\\!\\\\!\\\\!\\perp \\operatorname{PA}$. We will amend this in the final version.\n\n> ### **Q7: Intuitively, in section 6.1, given a point observation $X=x$, how possible could we learn the whole ellipse? Could you elaborate more on the experiment and connect it to the theoretical assumptions and guarantees?**\n\nA point observation of $X = x$ fixes the location but underdetermines the shape of the ellipse; thus, infinitely many ellipses can be drawn through that point. By further observing the angle $\\operatorname{PA} = \\operatorname{pa}$, it is possible to infer the semi-major and -minor parameters of the ellipse via causal *abduction*. The whole ellipse for an observed $(\\operatorname{pa}, x)$ can then be estimated by answering counterfactual queries for all angles in $(0,2\\pi)$.\n\nWe can provably learn the correct model by using our proposed Markovian OT-Flow coupling. In contrast, the Naive coupling violates Markovianity, entangling $U$ with $\\operatorname{PA}$, and yielding inconsistent counterfactuals. This constructed scenario is crucial for verifying our theoretical claims. It shows that when our assumptions and constraints hold (i.e. Markovian/FC setting, monotone OT-Flow map), counterfactual identification in $d>1$ from observational data is possible.\n\nMehta, R. CF-Seg: Counterfactuals meet Segmentation. arXiv preprint arXiv:2506.16213."}, {"review_id": "igRo93nCEk", "reviewer": "Reviewer_5cpq", "summary": "This work provides theoretical guarantees for counterfactual identifiability in structural causal models (SCMs). The main result of the paper is an extension of Theorem 5.1 of Nasr-Esfahany et al. (2023) to the multivariate setting. The authors also provide empirical evaluations of their methodology.\n\nMain assumptions: \n* A SCM with independent noise variables, that is, for all $i \\in \\{1, \\dots, n\\}$ \n\t$$X_i = f_i(X_{pa(i)}, U_i) \\qquad U_i \\perp X_{pa(i)}$$\nwhere $U_i$ and $X_i$ are assumed to be $d$-dimensional random vectors.\n* For all $i \\in \\{1, \\dots, n\\}$ the functions $f_i$ are bijective with respect to $U_i$.\n* Knowledge of the observational distribution of the SCM and its induced graph is assumed.\n\nGoal:\nThe goal is counterfactual identifiability, that is, the authors want to provide conditions under which a counterfactual query $Q$ is identifiable. \n\nUnder the aforementioned assumptions and in the univariate case, that is $d=1$, and assuming strict monotonicity of $f_i$ Nasr-Esfahany et al. (2023) establish counterfactual identifiability. The goal of the present paper is to extend this result to the multivariate case, that is, $d>1$.", "questions": "* Do you require strong monotonicity of the $f_i$'s? If not, how are your assumptions distinguished from strong monotonicity?\n\n* The proof of Proposition A.8 is not quite clear to me. In the proposition, it says that there is a function $g$ for every combination of $u$ and $pa$. In general these functions can be different. It is therefore unclear to me how the authors arrive at equation (60) as $g_{u_2}$ and $g^{-1}_{u_1}$ seem to not necessarily be equal.\n\n* If I understand Theorem 5.1 of Nasr-Esfahany et al. (2023) correctly, their results extends to the multivariate case when there is one noise variable for each $X_i$. In this context, it would be interesting to know what the significance is of modelling multiple noise variables for a single observed variable? Are there specific applications you have in mind? For example, intuitively, why would multiple noise variables be advantageous for the X-Ray experiment?", "limitations": "It would be nice to have a discussion about what it means that $X_i$ and $U_i$ are assumed to have the same dimension $d$, and about the implications of the bounded support assumption for the noise $U_i$.", "rating": 3, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 2, "significance": 2, "originality": 3, "strengths_and_weaknesses": "Quality:\nI am not familiar enough with the optimal transport literature to verify that the proofs are correct. However, it is sometimes not clear to me what assumptions are made to get the stated conclusion. Specifically, for Proposition 4.5 and 4.11. For example, it is not clear to me what “satisfying a proposition” (line 181, 228) means. I am also unsure about the correctness of Proposition A.8. I go into more detail about the proof of Proposition A.8 in the section \"Questions\".\nI like the empirical evaluations; they are clear and seem practically relevant.\n\nClarity:\nIn my opinion clarity is a significant issue with this paper. Especially in the theoretical part, the paper is missing mathematical clarity.\nIt is often unclear what assumptions are made. For example, it is unclear to me what precise assumptions are made in Proposition 4.5 and Proposition 4.11.  Do we assume monotonicity of $f_i$? It would also be good to see a discussion about the assumptions and what they mean intuitively. Below is a list of sections where I believe there are issues with clarity:\n* The term \"counterfactual query\" is used in Definition 4.2 but never defined mathematically. \n* Please only number the equations that you are actually referencing.\n* It seems like the counterfactual dynamic OT map $T^* $ is used in the main text but not defined. Maybe Equation (77) in the Appendix defines $T^*$? If so, please state this explicitly in the main text.\n* Sometimes \"$\\coloneqq$\" is used for definitions and sometimes not. Please be consistent in your use of notation.\n* Definition 4.2 is a definition and a statement in one. Please use the definition environment for definitions only.\n* $T$ seems to be a function $T: \\mathbb{R}^{d} \\times \\mathbb{R}^k \\to \\mathbb{R}^d$ (line 187) and $T: \\mathbb{R}^{d} \\times \\mathbb{R}^{d} \\times \\mathbb{R}^k \\to \\mathbb{R}^d$ (for example line 181). If you abuse notation, please indicate this clearly.\n* It also seems like $T$ is defined outside an environment and inside. Specifically, in line 187, 193 and 227. Are these all the same function?\n* There is no such thing as “satisfying a proposition” (line 181, 228 for example). Please state clearly what assumption you make.\n* Sometimes you use $U$ and $X$ and sometimes $\\boldsymbol{U}$ and $\\boldsymbol{X}$. It seems like $U$ and $X$ were never defined. It is unclear to me whether both notations mean the same thing.\n* There is a “motivating example” in the Appendix under “Proofs: Counterfactual Identifiability”. This seems misplaced.\n* In Equation (9) on the left hand side $x$ seems to be fixed while on the right hand side this is not the case.\n* “… absolutely continuous with densities strictly positive on a bounded, open, convex domain w.r.t. the Lebesgue measure”. It probably should be “… absolutely continuous w.r.t. the Lebesgue measure and with densities strictly positive on a bounded, open, convex domain.”\n* Propositions in the Appendix do not have the same numbers as the corresponding Propositions in the main text.\n* In line 612 the authors require “smooth respective densities” which is an assumption that does not appear in the propositions of the main text.\n* Line 726: Does the fact that the “velocity field is Lipschitz continuous” follow from a previously made assumption? Because it sounds like it is something that is assumed on the spot. If so please include all such assumptions in the proposition.\n\nSignificance:\nI would have liked to see more discussion about the relevance of multivariate noise vectors for counterfactual identifiability. It seems to me that Theorem 5.1 of Nasr-Esfahany et al. (2023) trivially extends to the multivariate case if there is one noise variable for each $X_i$, that is, $X_i = f_i(X_{pa(i)}, U_i)$. For which problems is it important to model multiple noise variables per observed variable?\n\nOriginality:\nThis is difficult for me to answer as I do not fully understand the assumption made in this paper. However, it seems that the proof of Theorem 5.1 of Nasr-Esfahany et al. (2023) can be extended to the multivariate setting by assuming strong monotonicity of the $f_i$ (i.e. the generalisation of monotonicity to the multivariate setting) and then use Brenier’s Theorem to get the unique, monotone transport map that pushes $U$ to $X$. $T$ is then the replacement of the inverse of the CDF in the proof of Theorem 5.1 of Nasr-Esfahany et al. (2023). The theoretical results seem to mirror this idea. However, I cannot find the strong monotonicity assumption for the $f_i$. Could the authors clarify this point?", "was_revised": true, "final_justification": "I would like to thank the authors for the productive discussion.\n\nI acknowledge their commitment to revising the manuscript. However, the scope of the necessary changes seem significant, and I therefore believe the paper is not yet ready for publication without addressing two key areas.\n\n$\\textbf{1. Framing of the Contribution: }$ In my opinion the paper's central claim needs to be framed more accurately. My understanding is that Theorem 5.1 of Nasr-Esfahany et al. (2023) already handles the multidimensional case. The key contribution of this work, therefore, appears to be the generalization from using a single noise variable for each $X_i$ to using a multi-dimensional noise vector. This is an important distinction, and the justification for this specific extension is not yet sufficiently rigorous, as discussed in my previous point (Q7). The authors should clarify this point and ensure their claims precisely reflect their contribution.\n\n$\\textbf{2. Mathematical Clarity: }$ Regarding mathematical clarity, I appreciate the authors' commitment to addressing the points I raised. However, I believe the required changes are substantial and go beyond minor edits. The presentation currently lacks precision in several key places, which hinders a reader's ability to verify the technical contributions. Because significant parts of the argument may need to be rewritten for clarity, I think a thorough revision cycle would be more appropriate to ensure the paper is clear and verifiable.\n\nGiven the authors' constructive engagement, I will raise my score to 3.", "rebuttal": "We sincerely appreciate the reviewer’s time in preparing their feedback. We are pleased that the empirical evaluations were found to be **\"clear and practically relevant**\". Although we respectfully disagree that clarity is a significant issue, as corroborated by the other reviewers, we commit to incorporating all the suggested changes to further improve clarity. We consider these to be minor and straightforward to address. Please find below our detailed responses to each concern.\n\n> ### **Q1: It is unclear to me what precise assumptions are made in Proposition 4.5 and Proposition 4.11. Do we assume monotonicity of f_i? It would also be good to see a discussion about the assumptions and what they mean intuitively.**\n\nFor an intuitive discussion about the assumptions, we invite the reviewer to refer to the dedicated *Assumptions & Plausibility* section in **Appendix A.5**, which appears to have been missed.\n\nTo clarify, we treat 'vector-monotonicity' (c.f. **Eq. (36)**) as a maintained hypothesis about the true data-generating function when applying the model to data. However, it was not obvious that the counterfactual transport map $T^*$ would be strictly vector-monotone in $x$, and we do not assume it a priori in the theory; we prove it is true given only the standard assumptions outlined below. \n\nThe assumptions made to prove **Proposition 4.5** are stated at the beginning of the proposition, precisely:\n1. Equal dimension $\\operatorname{dim}(X) = \\operatorname{dim}(U)$;\n2. Markovian setting $U \\perp\\\\!\\\\!\\\\!\\perp PA$;\n3. $P_U$ and $P_{X|\\mathbf{PA}}$ are absolutely continuous w.r.t. the Lebesgue measure, with strictly positive and bounded densities on a bounded, open, convex domain.\n\nItem 3 above corresponds to the assumptions needed for Caffarelli (1992)’s regularity theorem, so the result applies directly, ensuring that the OT map is bijective. For a helpful summary of Caffarelli (1992)’s results, the reviewer may refer to Theorem 12.50 in Villani (2008). In combination, our assumptions ensure that $P_U$ and $P_{X|\\mathbf{PA}}$ admit well-behaved densities and guarantee the existence and uniqueness of a dynamic OT map (Brenier 1991; Caffarelli 1992).\n\nThe assumptions made to prove **Proposition 4.11** are stated at the beginning of the proposition, namely:\n\n1. Assumptions 1, 2 and 3 as listed above for Proposition 4.5;\n2. $P_U$ is the continuous uniform measure on $[0,1]^d$;\n3. $T$ is the time-1 dynamic OT map that pushes $P_U$ forward to $P_{X|\\mathbf{PA}}$ described in Proposition 4.5.\n\nWe then prove that the counterfactual transport map $T^*$ must be strictly vector-monotone in $x$. To the best of our knowledge, this is a new result, and it extends the classical monotone-quantile notion to multivariate settings ($d > 1$), thereby guaranteeing counterfactual identifiability from observational data. \n\n> ### **Q2: $T$ seems to be a function $T : \\mathbb{R}^d \\times \\mathbb{R}^k \\to \\mathbb{R}^d$ (line 187) and $T : \\mathbb{R}^d \\times \\mathbb{R}^d \\times \\mathbb{R}^k \\to \\mathbb{R}^d$ (for example line 181). If you abuse notation, please indicate this clearly.**\n\nThank you for raising this. Firstly, we presume the reviewer meant to state that the map on **line 181** is $T : \\mathbb{R}^k \\times \\mathbb{R}^k \\times \\mathbb{R}^d \\to \\mathbb{R}^d$  rather than $T : \\mathbb{R}^d \\times \\mathbb{R}^d \\times \\mathbb{R}^k \\to \\mathbb{R}^d$. To clarify, we intended to redefine $T$ on **line 187** but recognise that this can be made clearer by avoiding an overload in notation with the map on **line 181**, for example. \n\nTo resolve this ambiguity, we will define the counterfactual transport map using $T^* : \\mathbb{R}^k \\times \\mathbb{R}^k \\times \\mathbb{R}^d \\to \\mathbb{R}^d$ before its first use in **Eq. (8)**, also ensuring consistency with its use in **Proposition 4.11**.\n\n> ### **Q3: Sometimes you use $X$ and $U$ and sometimes $\\boldsymbol{X}$ and $\\boldsymbol{U}$. It seems like $X$ and $U$ were never defined. It is unclear to me whether both notations mean the same thing.**\n\nAs per **Definition 3.1**, we use $\\mathbf{X}$ and $\\mathbf{U}$ to denote the set of all the endogenous and exogenous variables in the SCM, respectively. Without loss of generality, we then define $X$ and $U$ on line 142, as relating only to the $i^{\\text{th}}$ causal mechanism in the SCM. \n\n> ### **Q4: In line 612 the authors require “smooth respective densities” which is an assumption that does not appear in the propositions of the main text.**\n\nThank you for spotting this. We inadvertently stated this in the proof, but it is not necessary to obtain this result, so we will omit it in the final version. Please refer to our response to query **Q1** above for an intuitive discussion of the assumptions required to obtain the result.\n\n\n> ### **Q5: Does the fact that the “velocity field is Lipschitz continuous” follow from a previously made assumption? Because it sounds like it is something that is assumed on the spot. If so please include all such assumptions in the proposition.**\n\nThank you for the suggestion. We had treated Lipschitz continuity as a standard standing assumption in the Neural ODE literature. As in R. T. Chen et al. (2018) and Kidger et al. (2020), this condition holds when the velocity field is parameterised by a neural network with finite weights and Lipschitz activations (e.g. ReLU, Leaky-ReLU, tanh, softplus etc). Other smooth activations with bounded derivatives, such as GELU and SiLU, have finite global Lipschitz constants and can be scaled to 1-Lipschitz without changing qualitative behaviour. \n\nTo avoid ambiguity, we will state the above explicitly and explain that it follows immediately from our network parameterisation.\n\n> ### **Q6: The proof of Proposition A.8 is not quite clear to me. In the proposition, it says that there is a function g for every combination of u and pa. In general these functions can be different. It is therefore unclear to me how the authors arrive at equation (60) as $g_{u_2}$ and $g_{u_1}^{-1}$ seem to not necessarily be equal.**\n\nThank you, we are happy to clarify. Our **Proposition A.8** result mirrors Nasr-Esfahany et al. (2023)’s Proposition 6.2, but here it is extended to an optimal transport context.\n\nFirst, we prove that two equivalent transport maps $T^{(1)} \\sim_{\\mathcal{L}_3} T^{(2)}$, in the sense of **Definition 4.7**, produce the same counterfactuals. Transport equivalence implies that there exists an invertible transition map $g$ between the exogenous variables of the two transport maps $T^{(1)}, T^{(2)}$. That is, we have that:\n\n$u^{(1)} = g(u^{(2)}), \\quad \\text{and} \\quad u^{(2)} = g^{-1}(u^{(1)})$.\n\nThen, by the equivalence condition in **Definition 4.7**, the equation in question (i.e. **Eq. (60)**) follows directly by substitution:\n\n$T^{(1)}(g(u^{(2)}); \\mathbf{pa}^\\*) = T^{(2)}(g^{-1}(u^{(1)}); \\mathbf{pa}^\\*)$.\n\n\n> ### **Q7: If I understand Theorem 5.1 of Nasr-Esfahany et al. (2023) correctly, their results extends to the multivariate case when there is one noise variable for each x_i. In this context, it would be interesting to know what the significance is of modelling multiple noise variables for a single observed variable?**\n\nAs pointed out by Nasr-Esfahany et al. (2023) in their Section 5.1:\n\n*“It is not clear how to generalize the monotonicity condition to BGMs with multidimensional $V$, which is a known issue in Markovain causal structures (Nasr-Esfahany & Kiciman, 2023).”*\n\nOur contributions include precisely this missing generalisation to $d > 1$ for the Markovian, Instrumental Variable and Frontdoor Criterion settings. This unlocks counterfactual identifiability for multi-dimensional variables, such as images. \n\nThe Markovian setting in particular was an important open problem with close connections to both the (non-linear)-ICA (Hyvärinen and Pajunen 1999) and disentanglement (Locatello et al. 2019) literature, rendering our theoretical results broadly relevant and widely applicable. In the context of the most closely related work by Nasr-Esfahany et al. (2023), we can argue that our contribution represents a breakthrough in causal modelling, with significant implications for practical applications involving counterfactual image generation.\n\nThere seems to be a misunderstanding, as we never use multiple noise variables per observed variable. In all cases, we make the equal dimension assumption that $\\operatorname{dim}(X) = \\operatorname{dim}(U)$, for any multi-dimensional variable $X$ within an SCM $\\mathfrak{C}$ (see line 142). For an intuitive discussion about this assumption, please refer to the *Assumptions & Plausibility* section in **Appendix A.5**, and the response below.\n\n> ### **Q8: It would be nice to have a discussion about what it means that $X_i$ and $U_i$ are assumed to have the same dimension $d$, and about the implications of the bounded support assumption for the noise $U_i$.**\n\nThe equal-dimension assumption $\\operatorname{dim}(X) = \\operatorname{dim}(U)$ is commonplace in flow-based generative mode, as it ensures that a bijective transport map is mathematically possible. When this assumption does not hold for the true data-generating process, a bijective model can still be used by, for instance, augmenting the lower-dimensional variable with dummy coordinates. \n\nLastly, $P_U$ must be absolutely continuous w.r.t. the Lebesgue measure, with a strictly positive and bounded density on a bounded, open, convex domain for Caffarelli (1992)’s regularity theorem to apply. The implication of this is that the OT map is bijective provided that the target measure $P_{X|\\mathbf{PA}}$ meets the same criteria. We will make this clearer in the final version.\n\n**References**\n\nVillani, C., 2008. Optimal transport: old and new (Vol. 338, pp. 129-131). Berlin: springer.\n\nKidger, P., Morrill, J., Foster, J. and Lyons, T., 2020. Neural controlled differential equations for irregular time series. Advances in neural information processing systems, 33, pp.6696-6707."}, {"review_id": "bnJiWw1Krq", "reviewer": "Reviewer_LANz", "summary": "The work studies counterfactual identifiability from observational data in models with no unbserved confounders and high-dimensional variables. The proposed approach is is a dynamic optimal transport flow mechanism strictly monotone in the target variable. On the practical side, the authors describe a specialized Batch OT coupling approach suitable for Markovian model, i.e., exogeneity is satisfied. Extensive experiments are performed on real data to validate the the paper's claims.", "questions": "Please see weaknesses from above.\n\nAdditionally,\n\n(Q1) Can the results from this work be easily extended to more complex counterfactuals, e.g., nested counterfactuals?\n(Q2) Can you comment on the potential connection/overlap of the result from Proposition 4.3 and counterfactual identifiability results under *noise monotonicity* (Lu et.al. Triantafyllou et.al.).\n\nLu, Chaochao, et al. \"Sample-efficient reinforcement learning via counterfactual-based data augmentation.\" arXiv preprint arXiv:2012.09092 (2020).\nTriantafyllou, Stelios, et al. \"Agent-specific effects: A causal effect propagation analysis in multi-agent MDPs.\" ICML, 2024", "limitations": "Mentioning more clearly assumptions in conclusion and how they limit the pracictability of the work would be beneficial for readers.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 4, "significance": 3, "originality": 2, "strengths_and_weaknesses": "Strengths\n\n- The paper combines ideas from several fields and provides a rigorous theoretical analysis. I did not manage to follow all ideas and results presented, but regardless I recognize that this paper is well-written and nicely structured. Necessary background is well-covered given the space constraints. The authors have also made a very meticulous job in covering related work.\n- Experimental results are the main strength of this paper. Very cleverly executed and convincing in validating theoretical claims.\n\nWeaknesses\n\n- (Contributions) It feels that most results are derived or similar to prior work. I am familiar with the work from Nasr-Esfahany et al. (2023) which I understand does similar work for models with one-dimensional variables. What are the challenges going from the results from this work to the multi-dimensional problem? Can the authors highlight which contributions are novel and to what extent?\n- (Assumptions) To my understanding, the proposed approach does the following assumptions: exogeneity, dim(endogenous) = dim(exogenous), and exogenous distribution is uniform. It seems that this setting is too restrictive. Could the authors comment on this?", "was_revised": true, "final_justification": "Concerns about novelty were relieved.", "rebuttal": "We thank the reviewer for their positive assessment, noting our integration of ideas across fields, rigorous theoretical claims, clear structure, and thorough coverage of related work. We also appreciate their view that the experiments are a core strength: \"**very cleverly executed and convincing in validating theoretical claims**\". Please find our responses to all the reviewer's comments below.\n\n> ### **Q1: I am familiar with the work from Nasr-Esfahany et al. (2023) which I understand does similar work for models with one-dimensional variables. What are the challenges going from the results from this work to the multi-dimensional problem? Can the authors highlight which contributions are novel and to what extent?**\n\nWe thank the reviewer for their feedback and for the opportunity to clarify our contributions.\n\nAs pointed out by Nasr-Esfahany et al. (2023) in their **Section 5.1**:\n\n*“It is not clear how to generalize the monotonicity condition to BGMs with multidimensional $V$, which is a known issue in Markovian causal structures (Nasr-Esfahany & Kiciman, 2023).”*\n\nThis generalisation to $d > 1$ in Markovian causal structures was indeed an open problem, with close connections to both the (non-linear)-ICA (Hyvärinen and Pajunen 1999) and disentanglement (Locatello et al. 2019) literature. In multi-dimensional settings, one’s model is subject to indeterminacies, partly caused by symmetries of the prior $P_U$, which are unresolvable without further constraints.\n\nTo the best of our knowledge, we are the first to identify the constraints needed to provably solve this problem. We develop a new theoretical foundation for counterfactual identification of multi-dimensional variables by connecting dynamic Optimal Transport (OT) and causality. Building on this foundation, we introduced a  Markovian OT coupling flow which enables the estimation of high-dimensional, unique and monotone causal mechanisms. We demonstrate that our method enables consistent counterfactual inferences of high-dimensional variables ($d > 1$), given only observational data.\n\nEstablishing counterfactual identification of high-dimensional Markovian SCMs is both timely and important for the responsible use of causal AI, as many recent works rely on such assumptions to support causal claims (Pawlowski et al. 2020; Sanchez and Tsaftaris 2021; Sánchez-Martin et al. 2022; Monteiro et al. 2023; Ribeiro et al. 2023; Y. Chen, L. Tian, et al. 2024).\n\nAlthough not the main focus of our analysis, in **Appendix A.4** we also extend counterfactual identifiability results from observational data to non-Markovian multi-dimensional ($d > 1$) settings when the following common criteria apply: **(i)** Instrumental Variable (IV), **(ii)** Backdoor Criterion (BC), and **(iii)** Frontdoor Criterion (FC). We prove counterfactual identifiability when the FC applies under similar conditions to the BC setup in Nasr-Esfahany et al. (2023), and expand the instrumental variable results to $d > 1$. These results show our method is widely applicable and useful for all these common causal structures.\n\nTo conclude, in the context of the most closely related work by Nasr-Esfahany et al. (2023), we can argue that our contribution represents a breakthrough in causal modelling, with significant implications for practical applications involving counterfactual image generation.\n\n> ### **Q2: To my understanding, the proposed approach does the following assumptions: exogeneity, dim(endogenous) = dim(exogenous), and exogenous distribution is uniform. It seems that this setting is too restrictive. Could the authors comment on this?**\n\nWe thank the reviewer for raising this point, and invite them to refer to the dedicated section on *Assumptions & Plausibility* in **Appendix A.5**. In addition, we provide a more succinct response below.\n\nWe would like to point out that, aside from the uniform prior, these assumptions are commonplace in the literature (Pawlowski et al. 2020; Monteiro et al. 2023; Ribeiro et al. 2023; Nasr-Esfahany et al. 2023). What sets our work apart is that we provide counterfactual identifiability guarantees for the multi-dimensional case, a set of results that was previously missing for the Markovian, Instrumental Variable and Frontdoor criterion cases. These results are of vital importance to ensure that any causal claims derived from our models are properly supported.\n\nWe only require that the prior $P_U$ be a uniform measure for our most strict result. Under different priors, a counterfactual equivalence relation is obtained, in the spirit of many existing identifiability results (Khemakhem et al. 2021; Nasr-Esfahany et al. 2023; Hyvärinen et al. 2024). Thus, relaxing the uniformity constraint on the prior is permissible in exchange for less strict guarantees.\n\nWe close with Pearl and Bareinboim (2022):\n\n*\"Assumptions are self‑destructive in their honesty. The more explicit the assumption, the more\ncriticism it invites. [...] Researchers therefore prefer to declare 'threats' in public and make assumptions in private.\"*\n\nOur work counters this practice: we intend to make our assumptions and constraints explicit so they can be better understood, challenged, and ultimately relaxed in light of new evidence or insight.\n\n> ### **Q3: Can the results from this work be easily extended to more complex counterfactuals, e.g., nested counterfactuals?**\n\nIn **Appendix A.4.3**, we prove counterfactual identifiability for ($d>1$)-dimensional variables when the Fontdoor Criterion (FC) applies. In the FC setting, we already compute nested or cross-world counterfactuals of the sort found in mediation analysis, as described below.\n\nFor example, consider the graph: $X \\to M \\to Y$, with unobserved confounder $Z$, such that:  $X \\leftarrow Z \\to Y$.  \nTo compute a cross-world counterfactual, we nest two operations.\n\nGiven observed evidence $\\\\{M = m, X = x\\\\}$, we compute the counterfactual mediator $M_{x'}$ had $X$ been $x'$:\n\n$M_{x'} \\mid \\\\{M = m, X = x \\\\}$.\n\nIn structural equation form, using abduction-action-prediction, we compute:\n\n$m' = f_M(x', f_M^{-1}(x, m))$.\n\nWe then compute the counterfactual outcome (nested) given the counterfactual mediator:\n\n$Y_{x',M_{x'}} \\mid \\\\{ Y = y, M = m, X = x \\\\} $,\n\nwhich in functional form is given by the nested operation:\n\n$y' = f_Y(m', f_Y^{-1}(m, y)) = f_Y(f_M(x', f_M^{-1}(x, m)), f_Y^{-1}(m, y))$\n\nSo, yes, our approach already works for complex counterfactuals.\n\n> ### **Q4: Can you comment on the potential connection/overlap of the result from Proposition 4.3 and counterfactual identifiability results under noise monotonicity (Lu et.al. Triantafyllou et.al.).**\n\nWe thank the reviewer for the pointers. We will add a reference to Triantafyllou et al. (2024).\n\nDefinition 4.2 for noise monotonicity in Triantafyllou et al. (2024), as first specified by Pearl (2009) for binary SCMs, specifies a total order $\\leq$ on the domain of $V_i$, which naturally applies to scalars ($d=1$), but not vectors ($d>1$). In the absence of a meaningful natural order, readers have thus far only considered $V_i$ as a scalar quantity (Nasr-Esfahany et al. 2023), to the best of our knowledge.\n\nFor the $d=1$ case in **Proposition 4.3** (c.f. **Appendix A.2**), our result is practically equivalent to noise monotonicity, in line with Theorem 5.1 in Nasr-Esfahany et al. (2023). \n\nFor the $d>1$ case in **Proposition 4.3** (c.f. **Appendix A.2**), we provide a new result by identifying a generalised notion of monotonicity for vectors and its relation to Optimal Transport (OT). We find that since the OT map is the gradient of a convex function, it is a monotone operator, and thus the counterfactual transport map must be ‘vector-monotone’ in $x$; a desired property of counterfactual functions. This is not our main result, but it represents an important motivating step for our study of OT maps as causal mechanisms.\n\n**References**\n\nTriantafyllou, Stelios, et al. \"Agent-specific effects: A causal effect propagation analysis in multi-agent MDPs.\" ICML, 2024\n\nLu, Chaochao, et al. \"Sample-efficient reinforcement learning via counterfactual-based data augmentation.\" arXiv preprint arXiv:2012.09092 (2020).\n\nPearl, J. and Bareinboim, E., 2022. External validity: From do-calculus to transportability across populations. In Probabilistic and causal inference: The works of Judea Pearl (pp. 451-482)."}, {"review_id": "vwWEiYSRh3", "reviewer": "Reviewer_ugrj", "summary": "The authors consider counterfactuals of an outcome with respect to its parent set, assumed causally sufficient. Under strict monotonicity and regularity conditions referenced in the manuscript, the authors prove nonparametric identifiability from observational data using dynamic optimal transport. They construct an estimator and show that it provides consistent results in experiments (simulated and real world) where nonparametric baselines without rigorous identifiability proofs appear to be inconsistent.", "questions": "The authors' results apply to counterfactuals regarding atomic interventions and atomic conditions on evidences. Do the authors expect similar results to apply under soft interventions or distributional constraints on evidences?", "limitations": "Yes.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 4, "originality": 4, "strengths_and_weaknesses": "To the best of my knowledge, the authors are correct in saying that nonparametric counterfactual identification under monotonicity with multidimensional treatments has not been proven elsewhere. This result is therefore original and significant. The authors communicate their approach very well. \n\nThere are a few points which were not entirely clear to me: \n- The authors should expound upon the baselines and other estimators without consistency proofs: are these clearly incompatible the approach proposed by the authors?\n- While it is clear strict monotonically makes their approach easier to follow, where is it shown that bijectivity is insufficient alone? \n\nI disagree with the authors in the second paragraph of their introduction that point-identifiability is necessary for valid causal inference. Bounding, hypothesis testing and structure discovery are all valuable approaches in causal inference, even when a causal model is only partially identified. The authors make clear, however, that consistent counterfactual estimation requires point identification.", "was_revised": true, "final_justification": "Reviews have been largely positive with concerns raised by Reviewer LANz. Primarily, the reviewer is concerned as to the significance of an extension to multivariate outcomes. In my opinion this is very much a significant extension.\n\nEven forgoing counterfactual analysis and focusing on population-level treatment effects, I understand multidimensional outcomes in the instrument variables setup to significantly affect identifiability---in linear models we lose identifiability of the average treatment affect even with two outcomes. To my mind, the jump from single variable to multivariable outcomes is as significant as the jump from discrete/bounded outcomes to continuous unbounded ones.\n\nFor these reasons, I will be retaining my score.", "rebuttal": "We appreciate the reviewer's positive assessment and for recognising that our main theoretical result is \"**original and significant**\" and \"**has not been proven elsewhere**\". We devoted substantial effort to presenting all our results clearly, and we are pleased that this was recognised. We address all the reviewer's comments below.\n\n> ### **Q1: The authors should expound upon the baselines and other estimators without consistency proofs: are these clearly incompatible [with] the approach proposed by the authors?**\n\nWe appreciate the opportunity to elaborate on this. Baselines and other estimators without consistency proofs cannot safely be relied upon to support causal claims. In other words, a lack of identifiability, full or partial, implies the impossibility of making a precise causal claim. This means that the answer can be completely unrelated to the ground truth, and providing such guarantees is what causal methods are designed to do. \n\nFor example, the $\\mathrm{Naive}$ OT coupling flow baseline can model the observational distribution well, but the inferred counterfactuals are inconsistent, as the model violates Markovianity (c.f. **Figure 1**). Nasr-Esfahany et al. (2023) reported that their spline flow-based baseline also **failed** in the Markovian case, only succeeding when the backdoor criterion was applicable. In contrast, we demonstrate that our method produces consistent counterfactual inferences in both the Markovian and Frontdoor criterion settings (c.f. **Table 1**). We’ll clarify this further in the final version.\n\n> ### **Q2: While it is clear strict monotonically makes their approach easier to follow, where is it shown that bijectivity is insufficient alone?**\n\nThank you for pointing this out. We will make this clearer in the final version. \n\nIn short, a vector-valued function can be bijective but not monotone. As a simple example, consider the negative identity function:\n\n$T: \\mathbb{R}^d \\times \\mathbb{R}^k \\to \\mathbb{R}^d, \\quad T(u;\\mathbf{pa}) = -u$,\n\nwhich is bijective in $u$, but not vector-monotone:\n\n$ \\langle T(u_1;\\mathbf{pa}) - T(u_2;\\mathbf{pa}), u_1 - u_2 \\rangle = (-u_1 + u_2)^\\top(u_1 - u_2) = -\\lVert u_1 - u_2 \\rVert ^2_2 < 0$,\n\nas it violates the condition \n\n$ \\langle T(u_1;\\mathbf{pa}) - T(u_2;\\mathbf{pa}), u_1 - u_2 \\rangle \\geq 0$.\n\nAnother example is the $T^{(3)}$ mechanism of **Eq. (29)** in **Appendix A.1**, which is bijective but not monotone; thus, its counterfactual estimates are shown not to be consistent.  In **Appendix A.1**, we provide a gentle example showing that monotonicity is necessary for consistent counterfactual inferences in the Markovian setting.\n\n> ### **Q3: Bounding, hypothesis testing, and structure discovery are all valuable approaches in causal inference, even when a causal model is only partially identified. The authors make clear, however, that consistent counterfactual estimation requires point identification.**\n\n\nWe certainly agree that partial identification methods, such as Manski-style bounds (Manski 2003), constitute valuable causal inference approaches. Our intention was not to discount such efforts, but rather to stress that identification, partial or otherwise, is essential for valid counterfactual inference and is a prerequisite for making precise causal claims from observational data. Unfortunately, this fact remains underappreciated in ongoing empirical work on high-dimensional counterfactual inference, posing an operational risk which we highlight and seek to address. \n\nTo avoid misunderstanding, we will add the following to the text:\n\n*“In settings where full identification is impossible, researchers can still draw valuable causal conclusions through bounding, hypothesis testing, or structure discovery. In this work, however, we concentrate on scenarios that require point-identification to obtain consistent counterfactual estimates from observational data alone.”*\n\n> ### **Q4: The authors' results apply to counterfactuals regarding atomic interventions and atomic conditions on evidences. Do the authors expect similar results to apply under soft interventions or distributional constraints on evidences?**\n\nThank you for the insightful question. While further theoretical work is needed to accommodate *arbitrary* soft interventions, as they can alter the underlying causal mechanism and potentially violate (vector) monotonicity, the current method supports a suitably well-defined subset of soft interventions so long as they preserve those properties. \n\nSimilarly, distributional constraints on evidence (e.g. $ X  \\sim Q(X) $) are more general than atomic observations $X = x$, and may not yield point-identifiable counterfactuals without additional assumptions. We see this as fertile ground for future research and expect optimal transport to play a valuable role in the analysis. We will highlight this further in our Discussion section.\n\n**References**\n\nManski, C.F., 2003. Partial identification of probability distributions. New York, NY: Springer New York."}]

Author rebuttal:

Meta-review: This paper provides the counterfactual identifiability results, generalising the result of Nasr-Esfahany et al. (2023) to the multi-dimensional setting, using dynamic Optimal Transport (OT) approach. The contribution is a non-trivial advancement of the counterfactual inference for high-dimensional data such as images.

A majority of reviews have been largely positive, although there is a disagreement in terms of significance, framing of the contribution, and mathematical clarity. Reviewer `5cpq` questioned the significance of the contribution, citing that "Theorem 5.1 of Nasr-Esfahany et al. (2023) already handles the multidimensional case." and argued that the scope of the necessary changes seems too significant for the paper to be ready for publication. On the other hand, Reviewers `hTdz`, `LANz`, `ugrj` seem to agree that the contribution is significant enough for the paper to be published. Reviewer ugrj has specifically mentioned that "In [his] opinion this is very much a significant extension [of Nasr-Esfahany et al. (2023)]." 

After reading the paper, I concur with the majority of the reviewers that the contribution of this paper is substantial enough to warrant its publication at NeurIPS as a poster. This decision is based on the argument that Theorem 5.1 of Nasr-Esfahany et al. can already handles the multidimensional case. Nevertheless, this doesn't significantly diminish the contribution of this paper.
Decision: Accept (poster)