1 Introduction↩︎

Entropic Optimal Transport (EOT) with quadratic cost is well-established mathematical framework with strong theoretical properties that has found wide application in generative modeling and especially for unpaired domain translation. Despite this, the practical adoption of entropic transport methods has been limited by the lack of efficient and flexible algorithms. Existing approaches typically suffer from at least one of the following drawbacks: they are not amenable to simulation-free training, requiring costly sampling at each iteration [1]; they require adversarial optimization [2]; they involve training a sequence of models rather than a single objective [3]; they impose restrictive parametric forms on the transport plan [4]; or they turn out to be too sensitive to the entropic regularization strength [5].

In our paper, we make a decisive step towards solving the drawbacks of existing EOT methods, and propose a novel simulation-free training solver based on an innovative variational reformulation of the weak dual EOT objective. We present the following main contributions:

1. Variational dual objective (3.1). We derive an equivalent reformulation of the weak dual objective of entropic OT with quadratic cost in which the intractable log-partition term is replaced by an exact variational minimization over an auxiliary log-normalizer.

2. Simulation-free training solver (3.2). Building on this reformulation, we propose a simulation-free training variational solver that jointly learns the dual potential and the auxiliary normalizer via neural parameterizations (no MCMC).

3. Learning guarantees (3.3). We provide finite-sample learning guarantees for recovery of the entropic OT plan, decomposing error into estimation and approximation terms, and show vanishing approximation error under universal function approximation.

4. Evaluation (5). We evaluate our solver on synthetic data and unpaired image-to-image translation. For the latter, we highlight gains by comparing against other solvers that optimize the same weak dual objective.

2 Background↩︎

This section provides the necessary background on entropic optimal transport and its optimization using weak dual reformulation. In 2.1, we recall the entropic optimal transport problem with quadratic cost, introduce its weak dual formulation, and describe the structure of the optimal transport plan induced by the optimal dual potential. In 2.2, we clarify our learning setup. In 2.3, we review existing solvers based on weak dual objective, highlighting two main paradigms: simulation-based optimization via implicit sampling proposed by the authors of EgNOT [1] and simulation-free methods based on restrictive parametric assumptions proposed by the authors of LightSB [4]. This discussion motivates the need for a dual formulation amenable to simulation-free training while remaining sufficiently expressive, which we address in subsequent 3 by proposing our novel VarEOT solver.

2.1 Entropic Optimal Transport with the quadratic cost↩︎

Let \(p_0, p_1 \in \mathcal{P}_{\mathrm{ac}}(\mathbb{R}^{D})\) be the set of absolutely continuous Borel probability measures on \(\mathbb{R}^{D}\). Let \(\Pi(p_0,p_1)\) denote the set of couplings (transport plans) on \(\mathbb{R}^D\times\mathbb{R}^D\) with marginals \(p_0\) and \(p_1\). We write \(\pi(x_0,x_1)\) for a plan density; \(H\) is the differential entropy. The Entropic Optimal Transport (EOT) problem is given by: \[\begin{align} \mathrm{EOT}_{\varepsilon}(p_0,p_1) &\overset{\mathrm{def}}{=} \min_{\pi\in\Pi(p_0,p_1)} \Bigg\{ \underbrace{\underset{\pi(x_0, x_1)}{\mathbb{E}}\!\left[\frac{\|x_0-x_1\|^2}{2}\right]}_{\text{optimal transport term}} \nonumber \\ &\underbrace{- \varepsilon\int_{\mathbb{R}^D} \!\!H \!\big(\pi(\cdot\mid x_0)\big)p_0(x_0)dx_0}_{\text{entropic regularization term}} \Bigg\}. \label{eq:eot} \end{align}\tag{1}\] The entropic regularization term in 1 is due to [1]. There are equivalent forms of the term [6], [7] which differ only by additive constants that do not affect the solution.

The EOT problem admits a unique minimizer \(\pi^\ast\), referred to as the EOT plan. While the primal formulation 1 is conceptually appealing, it is computationally inconvenient, since enforcing the marginal constraints \(\pi \in \Pi(p_0,p_1)\) requires optimizing over a complex set of probability measures that does not admit a straightforward parametrization.

Weak dual form of EOT. Objective 1 admits the following weak dual representation [1]: \[\sup_{f\in L_1(\mathbb{R}^D)} \bigg\{\underbrace{ {\underset{p_1(x_1)}{\mathbb{E}}}f(x_1) - \varepsilon {\underset{p_0(x_0)}{\mathbb{E}}}\log Z(f,x_0)}_{\stackrel{\text{def}}{=}\mathcal{L}(f)} \biggr\}, \!\! \label{eq:dual-final}\tag{2}\] where the \(\sup\) is taken over a function \(f\in L_1(p_1)\) satisfying \(f\equiv-\infty\) outside the \(\mathrm{supp}(p_1)\) and \[Z(f,x_0) \stackrel{\text{def}}{=} \int_{\mathbb{R}^{D}} \exp\!\left( \frac{f(x_1) - \tfrac{1}{2}\|x_0-x_1\|^{2}}{\varepsilon} \right) dx_1 \label{eq:partition}\tag{3}\] is the partition function.

Optimal transport plan. Let \(f^*\) be an optimizer of 2 , the corresponding optimal transport plan \(\pi^*\) [1] can be recovered from it. By the disintegration with respect to the source marginal \(p_0\) we have: \[\pi^*(x_0,x_1) = \pi^*(x_1 \mid x_0)\, p_0(x_0), \label{eq:joint95pi95star}\tag{4}\] Then the density of conditional distribution \(\pi^*(\cdot \mid x_0)\) in 4 is directly defined by \(f^*\): \[\pi^*(x_1 \!\!\mid\! x_0) \!=\! \frac{1}{Z(f^*,x_0)} \exp\!\left(\! \frac{f^*(x_1) \!-\! \tfrac{1}{2}\|x_0\!-\!x_1\|^{2}}{\varepsilon} \right). \label{eq:cond95pi95star}\tag{5}\]

2.2 Computational EOT setup↩︎

In practice, the source and target distributions, \(p_0, p_1\), as well as the EOT objective 1 , could be expressed and treated in different ways. To avoid possible misunderstanding, below we formalize our practical learning setup:

This setup falls within continuous OT, in contrast to discrete OT [6], [8], which is designed to compute one-to-one or one-to-many correspondence directly between the collections of provided source and target samples. As a result, discrete OT approaches do not naturally accommodate the out-of-sample estimation demanded by continuous OT. In our manuscript, we focus exclusively on continuous OT approaches, treating discrete OT as a considerably different direction.

2.3 Existing Weak Dual Formulation Solvers↩︎

The practical optimization of problem 2 remains challenging due to the presence of the partition function \(Z(f,x_0)\), which is, in general, intractable to compute exactly. Below we review two representative strategies to optimize the semi-dual objective: (i) general neural potentials with implicit sampling (EgNOT), and (ii) simulation-free objectives enabled by restrictive parametric transport families (LightSB).

EgNOT solver. The authors of EgNOT [1] solve 2 by parametrizing \(f_{\theta}\) by a neural network and deriving the gradient of the weak dual objective \(\mathcal{L}(f_{\theta})\): \[\begin{align} \nabla_{\theta}\,\mathcal{L}( f_\theta) = - {\underset{p_0(x_0)}{\mathbb{E}}}\left[{\underset{\pi_{\theta}(x_1|x_0)}{\mathbb{E}}}\left[ \nabla_{\theta} f_\theta(x_1)\right]\right] \, \\ + {\underset{p_1(x_1)}{\mathbb{E}}}\left[\nabla_{\theta} f_\theta(x_1)\, \right], \end{align}\]
where \(\pi_{\theta}(x_1|x_0)\) is given by: \[\pi_{\theta}(x_1 \!\!\mid\! x_0) \!=\! \frac{1}{Z(f_{\theta},x_0)} \exp\!\left(\! \frac{f_{\theta}(x_1) \!-\! \tfrac{1}{2}\|x_0\!-\!x_1\|^{2}}{\varepsilon} \right).\] While flexible, this approach is not simulation-free: each evaluation of the loss or its gradient requires sampling from the model distribution \(\pi_{\theta}(x_1 \mid x_0)\), which is itself defined implicitly through the neural potential \(f_{\theta}\). In practice, this sampling step is carried out using Markov chain Monte Carlo (MCMC) methods [9]–[11], such as Langevin dynamics. However, MCMC-based sampling can be computationally expensive, sensitive to hyperparameters, and slow to mix, especially in high-dimensional settings.

LightSB solver. An alternative strategy is proposed in LightSB [4]. The authors introduce adjusted potential \(v_{\theta}\) and parametrization: \[\pi_{\theta}(x_1|x_0)=\frac{\exp\big(\frac{\langle x_0,x_1\rangle}{\varepsilon}\big)v_{\theta}(x_1)}{c_{\theta}(x_0)}, \label{eq:lightsb-plan-parametric-full}\tag{6}\]
where, \(c_{\theta}(x_{0})\stackrel{\text{def}}{=}\int_{\mathbb{R}^{D}}\exp\big(\frac{\langle x_0,x_1\rangle}{\varepsilon}\big)v_{\theta}(x_1)dx_1\) is the normalization. They then consider the optimization problem:: \[\min_{v_{\theta}} \biggr\{ {\underset{p_0(x_0)}{\mathbb{E}}}\log c_{\theta}(x_0)-{\underset{p_1(x_1)}{\mathbb{E}}} \log v_{\theta}(x_1)\biggr\}. \label{eq:lightsb-main-obj-feasible}\tag{7}\]
This problem is equivalent to the problem 2 under a different parametrization, specifically: \[\begin{gather} v_{\theta}(x_1) = \exp(-\frac{\|x_1\|^{2}}{2\varepsilon})\exp(\frac{f_{\theta}(x_1)}{\varepsilon}); \nonumber \\ c_{\theta}(x_0) = \exp(\frac{\|x_0\|^{2}}{2\varepsilon})Z(f_{\theta}, x_0). \nonumber \end{gather}\]
To solve the problem in practice, the authors of LightSB circumvent the intractability of normalization \(c_{\theta}(x_0)\) (equivalent of \(Z(f,x_0)\), eq. 3 ) by directly parameterizing the potential \(v_{\theta}(x_1)\) in form of a Gaussian mixture: \[\begin{gather} v_{\theta}(x_1) =\sum_{k=1}^{K'} \alpha_k \,\mathcal{N}(x_1\vert r_k,\varepsilon S_k), \end{gather}\] where \(\theta\stackrel{\text{def}}{=}\{\alpha_{k},r_{k},S_{k}\}_{k=1}^{K'}\) are the parameters: \(\alpha_{k}\geq 0\), \(r_k\in\mathbb{R}^{D}\) and symmetric \(0\!\prec\! S_{k}\!\in\!\mathbb{R}^{D\times D}\). Such a parameterization restricts conditional density \(\pi(x_1 \!\mid\! x_0)\) to form: \[\begin{gather} \pi_{\theta}(x_{1}|x_0)=\frac{1}{c_{\theta}(x_0)}\sum_{k=1}^{K'}\widetilde{\alpha}_{k}(x_0)\mathcal{N}(x_{1}|r_{k}(x_0),\varepsilon S_{k}), \nonumber \end{gather}\] where \(\widetilde{\alpha}(x_0)\) and \(r_k(x_0)\) are some functions of \(\{\alpha_{k},r_{k},S_{k}\}_{k=1}^{K'}\) and \(c_{\theta}(x_0)\) is given in a closed form: \[\begin{gather} c_{\theta}(x_0) = \sum_{k=1}^{K'}\widetilde{\alpha}_k(x_0). \nonumber \end{gather}\] While this parameterization leads to a fully tractable and simulation-free objective, it restricts the expressiveness of the method, as it limits admissible transport plans to a rather narrow parametric family.

Summary. Existing weak dual solvers trade off between expressiveness and tractability: EgNOT supports flexible potentials but requires MCMC during training, while LightSB is simulation-free but restricts the conditional plan family. Below, we present our novel method, VarEOT, which takes the best of two worlds by enabling simulation-free training without restricting \(\pi(x_1\!\mid\! x_0)\) to a narrow parametric family.

3 Variational Entropic Optimal Transport↩︎

In this section, we introduce our variational approach to entropic optimal transport. In 3.1, we derive a new variational dual formulation of the EOT objective that replaces the intractable log-partition function with a tractable variational upper bound, yielding a fully differentiable and simulation-free training loss. Building on this formulation, 3.2 presents a practical learning algorithm based on neural parameterization of the dual potential and the auxiliary variational function, together with a Langevin sampling procedure for generating transport samples. We further position our method relative to existing dual solvers and highlight its practical advantages. In 3.3 we conduct the analysis of our method from the perspectives of statistical learning theory (finite sampling learning guarantees and approximation with neural networks). All proofs are provided in Appendix 7.

3.1 New Variational Dual Formulation of EOT↩︎

Our goal is to propose a weak dual solver that, unlike EgNOT and LightSB, does not require simulation during training and allows for expressive parameterization.

A key challenge in this setting is differentiating through the partition function \(\log Z(f,x_0)\): as a logarithm of an expectation, it cannot be unbiasedly estimated from finite samples in a straightforward way. To overcome this difficulty, we adopt a variational approximation for the logarithm, which allows us to construct a tractable, differentiable surrogate for the dual objective. The details of this variational approach are formalized in the proposition below.

Thanks to the obtained estimate, we can obtain tractable simulation-free loss:

This novel variational dual form allows us to overcome the original problem of estimation of log partition function \(\log Z(f, x_0)\) in the weak dual form 2 .

For convenience we define: \[\begin{align} \pi^{f,\xi}(x_0,x_1\!) &\!\stackrel{\text{def}}{=}\! \frac{(2\pi \varepsilon)^{-\frac{D}{2}} p_0(x_0\!)}{\exp(\xi(x_0))}\!\exp\!\bigg[\!\frac{f(x_1)\!-\!\frac{1}{2}\|x_0\!-\!x_1\|^2}{\varepsilon}\bigg]; \\ \pi^{f}(x_0, x_1) &\stackrel{\text{def}}{=}\pi^{f, \xi_f} (x_0, x_1), \end{align}\]
where \(\xi_f\) is due to ?? . The following theorem establishes that the VarEOT optimality gap directly corresponds to the KL discrepancy between the recovered measure \(\pi^{f,\xi}\) (\(\pi^f\)) and the EOT plan \(\pi^*\).

Our Theorem 3 certifies that optimizing objective 8 directly enables us to approximate the ground truth EOT plan. Additionally, eq. 9 suggests that at inference stage it is better to use \(\pi^f\), not \(\pi^{f, \xi}\), which we exploit in our practical implementation, see the next section 3.2.

3.2 Practical Algorithm↩︎

Training. The variational formulation derived in 2 leads to a fully tractable and simulation-free training optimization objective, given in 8 . In practice, we parametrize both the potential \(f : \mathbb{R}^D \to \mathbb{R}\) and the auxiliary variational function \(\xi : \mathbb{R}^D \to \mathbb{R}\) by neural networks, denoted by \(f_{\theta}\) and \(\xi_{\psi}\), respectively

The resulting training procedure consists of maximizing an empirical estimate of the variational dual loss with respect to both parameter sets \(\theta\) and \(\psi\). Crucially, both networks \(f_\theta\) and \(\xi_\psi\) are trained jointly in a single optimization loop, with no alternating or adversarial steps required. Given mini-batches of samples \(\{x_i^0\}_{i=1}^{N_0} \sim p_0\) from the source distribution and \(\{x_j^1\}_{j=1}^{N_1} \sim p_1\) from the target distribution, the expectations in eq. 8 are approximated using Monte Carlo sampling with i.i.d.Gaussian noise variables \(z_{i,k} \sim \mathcal{N}(0,I)\), yielding the empirical loss, up to an additive constant in 8 : \[\label{eq:empirical-loss} \begin{align} \widehat{\mathcal{L}}(f_{\theta}, \xi_{\psi}) \stackrel{\text{def}}{=} \frac{1}{N_1} \sum_{j=1}^{N_1}{f_{\theta}(x_j^1)} - \varepsilon\frac{1}{N_0} \sum_{i=1}^{N_0} \xi_{\psi}(x_i^0) - \\ \frac{\varepsilon}{N_0K} \sum_{i,k=1}^{N_0,K} \exp\!\left( \frac{f_{\theta}(x_i^0 + \sqrt{\varepsilon}\, z_{i,k})}{\varepsilon} - \xi_{\psi}(x_i^0) \right). \end{align}\tag{10}\] Importantly, this approximation does not require sampling from the model distribution itself, in contrast to energy-based approaches such as EgNOT. The complete training procedure is summarized in 1. To mitigate numerical instability arising from the exponential terms in 10 , we clip their values prior to gradient computation and additionally apply gradient clipping during training; see Appendix 8 for details.

Inference. Following the hint in equation 9 , after the training we only use the potential function \(\hat{f}_{\theta}\) for inference, whereas the function \(\hat{\xi}_{\psi}\) is not used. The corresponding entropic optimal transport plan is implicitly defined via the conditional distribution 5 . To generate samples from this conditional distribution, we employ Langevin dynamics targeting the unnormalized density \[\pi(x_1 \mid x_0) \propto \exp\!\left(\frac{ \hat{f}_{\theta}(x_1) - \frac{1}{2}\|x_1 - x_0\|^2}{\varepsilon} \right).\] The resulting sampling procedure is detailed in 2. In all experiments, we find that \(10^1\)-\(10^3\) Langevin steps with an appropriate step size suffice for high-quality samples.

To highlight the practical advantages of VarEOT, we compare it with two representative dual-solver approaches: EgNOT [1] and LightSB [4] (cf. 2.3). Table 1 summarizes their properties in terms of simulation-free training and structural constraints on the transport plan.

VarEOT combines the best of both worlds: it is simulation-free training, unlike EgNOT, and it imposes no restricted parameterization on the conditional plan, unlike LightSB.

3.3 Finite Sample Learning Guarantees↩︎

In this subsection we quantify the discrepancy between the transport plan recovered by VarEOT and the ground truth EOT solution. Our method does not have access to \((p_0,p_1)\) explicitly; instead it works with finite samples (see our computational setup §2.2) and with potentials restricted to parametric classes. This naturally introduces several sources of error in practice: finite-sample error (we only observe i.i.d. datasets of sizes \(N\) and \(M\) from \(p_0\) and \(p_1\)); function class restriction (we optimize over a class of neural networks rather than over all continuous potentials); optimization error (stochastic training, finite-time optimization). Our rigorous theoretical analysis below focuses on the first two items. The third item depends on the chosen optimizer and sampling scheme and is treated separately in practice.

Throughout, the function classes \(\mathcal{F}\) and \(\Xi\) are understood as classes of feedforward neural networks with biases, prescribed depth and width, uniformly bounded layer norms, and clipped outputs. Here, clipped outputs mean that the final network output is truncated to lie in a fixed bounded interval. We additionally assume that the hidden activation \(\sigma\) is \(L_\sigma\)-Lipschitz, non-polynomial and satisfies \(\sigma(0)=0\). Thus, \[\label{eq:classesFXI} \mathcal{F}\;:=\; \mathsf{NN}_\theta^{\sigma}(D,1), \quad \Xi \;:=\; \mathsf{NN}_\eta^{\sigma}(D,1),\tag{11}\] We begin with a decomposition of the excess \(\mathrm{KL}\) discrepancy between the ground-truth EOT optimizer \(\pi^\ast\) and the learned plan \(\pi^{\hat{f}}\) (cf. 9 ), where \[\begin{align} (\hat{f},\hat{\xi}) = \mathop{\mathrm{arg\,max}}_{f\in\mathcal{F},\;\xi\in\Xi}\widehat{\mathcal{L}}(f,\xi). \end{align}\] Here \(\xi\) enters only through the training objective, whereas the transport plan used at inference time is determined exclusively by the learned potential \(\hat{f}\); accordingly, the final discrepancy is measured in terms of \(\pi^{\hat{f}}\), cf. §3.2.

Proposition 4 separates the excess error into a statistical term (estimation error) and a representational term (approximation error). The first measures the uniform discrepancy between the empirical and population objectives over the prescribed neural-network classes, while the second quantifies the gap from restricting the optimization to these classes.

The next theorem controls the statistical term. It shows that, for fixed network classes, the estimation error decays at the usual sublinear rate as the number of samples increases. Importantly, decay rates do not depend on the dimension \(D\).

We next turn to the approximation term. The following statement asserts that it can be made arbitrarily small by choosing sufficiently expressive neural-network classes.

Thus, the bias introduced by the neural-network parameterization can be driven arbitrarily close to zero by enlarging the expressive power of the admissible classes.

Summary. Proposition 4, together with Theorems 5 and 6, yields the standard learnability picture for EOT plan recovery. For fixed neural-network classes, the estimation gap vanishes as the sample sizes grow. Meanwhile, by increasing the expressiveness of the classes, the approximation error can be made arbitrarily small. Consequently, with sufficiently many samples and sufficiently expressive neural-network parameterizations, the learned plan \(\pi^{\hat{f}}\) approaches the EOT solution in excess \(\mathrm{KL}\) discrepancy.

3.3.0.1 Relation to prior work.

Related error decompositions appear in [1], [14]. While [1] does not provide a fully detailed breakdown of the individual error terms for EgNOT, [14] develops a finer-grained statistical analysis but for a different objective (entropic barycenters). We follow a similar template, yet provide a detailed analysis for our more elaborate variationally derived VarEOT functional.

4 Related Work↩︎

In this section, we review the existing continuous EOT approaches, and briefly compare them with our solver. For clarity, we introduce a taxonomy of the existing approaches below, based on the utilized mathematical framework.

Weak dual EOT (a.k.a. semi-dual EOT) methods [1], [4] optimize the objective eq. 2 . They are the most similar to our VarEOT in terms of methodology. The detailed discussion of the approaches is given in 2.3. Crucially, our proposed method successfully overcomes the limitations of existing weak dual approaches while retaining their advantages.

Dual EOT [5], [15], [16] exploit an alternative dual form of the EOT. These methods simultaneously optimize a pair of dual potentials \((u, v)\) in a \(\max \max\) optimization procedure known as the Sinkhorn algorithm. While bearing certain resemblance to our VarEOT (e.g., simulation-free training), dual EOT approaches have pitfalls underscoring their performance:

• Dual potentials \((u, v)\) do not fully recover the (conditional) EOT \(\pi^*(\cdot \vert x_0)\). In particular, an auxiliary score-based model for the target distribution \(p_1\) is required to enable sampling from \(\pi^*(\cdot \vert x_0)\) [5]. Otherwise, there are only some heuristics (e.g., barycentric projection) that allow the recovered solution to be cast as a generative model capable of data-to-data translation [15], [16].

• Due to optimization peculiarities, dual EOT approaches tend to be unstable under small regularization coefficient \(\varepsilon\). This behavior is reported in [5] and supported by [1].

Overall, VarEOT provides a more user-friendly and production-ready framework with less engineering overhead required to adapt the method to applications.

Schrödinger bridge (SB) methods [2], [3], [12], [13], [17]–[21] cast the EOT in a dynamic manner and recover a solution in the form of a stochastic differential equation with learnable drift, with the endpoints given by source \(p_0\) and target \(p_1\) distributions. The underlying principles of SB-based approaches are different, e.g., adversarial optimization, gradual iterative markovian/proportional fitting procedure. However, compared to our method, majority of them are not simulation-free at training. The exceptions are [19] (based on mini-batch EOT approximation, which may not be accurate) and [21] (similar to [4] uses restricted Gaussian mixture approximation).

5 Experimental Illustrations↩︎

We evaluate VarEOT on both synthetic and real-world data. 5.1 presents illustrative two-dimensional experiments demonstrating the effect of the regularization parameter \(\varepsilon\) on the learned transport plan. 5.2 benchmarks VarEOT on high-dimensional biological single-cell data from the MSCI dataset, comparing against EOT/SB solvers. 5.3 comprehensively evaluates on unpaired image-to-image translation using the standard ALAE latent space protocol, comparing VarEOT against related entropic transport baselines across multiple translation tasks and regularization regimes. Technical details, including architecture specifications and hyperparameters, are provided in Appendix 8.

5.1 Two-dimensional Examples↩︎

We begin with a standard illustrative two-dimensional experiment that provides intuition into the role of the entropy regularization parameter \(\varepsilon\) in shaping the learned conditional transport plan \(\pi_{\theta}(x_1 \mid x_0)\). Specifically, we consider a setting in which a Gaussian source distribution is transported to a Swiss Roll target distribution. We apply our method for \(\varepsilon \in \{10^{-2},\,10^{-1},\,10^{0}\}\), and visualize the resulting transport maps in Fig. 3.

For small values of \(\varepsilon\), the learned transport is nearly deterministic, with little variation in sampling endpoints. As \(\varepsilon\) increases, the transport becomes progressively more stochastic: the endpoints diversify, and the conditional distributions \(\pi_{\theta}(x_1 \mid x_0)\) spread over a wider region of the target space.

5.2 Single Cell Data↩︎

One of the important applications of EOT is the analysis of biological single-cell data [22]–[24]. We evaluate VarEOT on the high-dimensional MSCI dataset (Kaggle competition “Open Problems — Multimodal Single-Cell Integration”) first used in [19]. The dataset consists of single-cell data from four human donors at 4 time points (days \(2\), \(3\), \(4\), and \(7\)). We solve the EOT problem between distribution pairs at days \(2\to4\) and \(3\to7\), and evaluate how well the solvers recover the intermediate distributions at days \(3\) and \(4\) respectively. We consider PCA projections with \(\text{DIM}=50, 100, 1000\) components and report energy distance [25] (ED) in Table 2. For VarEOT, the experimental setup (architecture, number of iterations, \(\varepsilon\), and learning rate) follows the EgNOT configuration from [4] with \(K{=}128\) particles; we refer the reader to Appendix D.4 of [4] for full details. VarEOT achieves quality comparable to Langevin-based and Minimax solvers for \(\text{DIM}=50\) and \(100\), while being \(5\)–\(15\times\) faster (\(3\) m vs.\(19\)–\(45\) m) on less powerful hardware (GTX 1080 vs.V100). For \(\text{DIM}=1000\), VarEOT matches the quality of the Langevin-based solver (\(1.47\) vs.\(1.46\)) with comparable training time (\(14.8\) m vs.\(15\) m). The details of preprocessing, hyperparameters and baselines are in Appendix 8.

5.3 Unpaired Image-to-Image Translation↩︎

Unpaired image-to-image translation is a standard benchmark in the entropic optimal transport (EOT) and Schrödinger Bridge (SB) literature [5], [26], [27], as it naturally combines high-dimensional data and the need for stochastic transport plans. In particular, translation in the latent space of pretrained autoencoders has become a popular and convenient setup for comparing entropic transport methods [1], [4], [28].

We follow the widely adopted ALAE protocol [4], [21], [28]–[31] based on ALAE autoencoder [32]. Specifically, we use a pretrained ALAE autoencoder trained on the full \(1024 \times 1024\) FFHQ dataset [33], which contains approximately \(70\)K human face images. The first \(60\)K images are used for training and are split into (male, female) and (child, adult) subsets. Using the fixed ALAE encoder, we extract \(512\)-dimensional latent representations \(\{z_0^n=\mathrm{Enc}(x_0^n) \}_{n=1}^N\) and \(\{z_1^m=\mathrm{Enc}(x_1^m)\}_{m=1}^M\) corresponding to the source and target domains. We consider \(4\) setups: male to female (M\(\to\)F), female to male (F\(\to\)M), adult to child (A\(\to\)C), and child to adult (C\(\to\)A).

Training. Given unpaired samples from the two latent distributions, we learn a latent entropic optimal transport plan \(\pi_\theta(z_1 \mid z_0)\) using our variational formulation. The model is trained entirely in latent space and does not require paired data or image-level supervision.

Inference. To translate a previously unseen image \(x_0^{\text{new}}\) (from the remaining \(10\)K test images), we (i) encode it as \(z_0^{\text{new}} = \mathrm{Enc}(x_0^{\text{new}})\), (ii) sample \(z_1 \sim \pi_\theta(z_1 \mid z_0^{\text{new}})\) using Langevin dynamics, and (iii) decode \(x_1 = \mathrm{Dec}(z_1)\).

Evaluation Protocol. We evaluate the quality of generated samples using two metrics: Fréchet Inception Distance (FID) [34] and Learned Perceptual Image Patch Similarity (LPIPS) [35]. FID is computed between the set of test translated images \(\{x_1^n = \mathrm{Dec}(z_1^n)\}_{n=1}^{N}\), where \(z_1^n \sim \pi_\theta(z_1 \mid z_0^n)\), and the set of ALAE-reconstructed target images \(\{\tilde{x}_1^m = \mathrm{Dec}(\mathrm{Enc}(x_1^m))\}_{m=1}^{M}\). LPIPS is measured between the source \(x_0\) and the translated \(x_1 = \mathrm{Dec}(z_1)\), where \(z_1 \sim \pi_\theta(z_1 \mid \mathrm{Enc}(x_0))\). It measures the input-output similarity.

5pt

Results. We comprehensively evaluate VarEOT across all four translation setups: M\(\to\)F, F\(\to\)M, A\(\to\)C, and C\(\to\)A. Figure 4 shows qualitative results for the M\(\to\)F setup with \(\varepsilon = 1.0\). We report translations obtained using different numbers of Langevin steps at inference time (NFE \(= 2, 10, 1000\)), demonstrating that VarEOT produces visually plausible and stable results even with a very small number of sampling steps. Additional results for alternative translation directions (F\(\to\)M, A\(\to\)C, C\(\to\)A) are provided in Appendix 9 (Figures 5, 6, and 7), confirming consistent behavior across all setups. Quantitative evaluation using FID and LPIPS metrics is summarized in Table 3, covering a wide range of regularization parameters \(\varepsilon \in \{0.1, 0.5, 1.0, 10.0\}\) and all four translation directions. We compare VarEOT against closely related entropic transport baselines: LightSB and EgNOT, showing competitive or improved performance across all considered regimes.

A qualitative analysis of the effect of the regularization parameter \(\varepsilon\) on sample diversity is provided in Figure 8, illustrating how varying \(\varepsilon\) controls the trade-off between transport cost and output stochasticity. Finally, Figure 9 provides an in-depth study of FID performance as a function of the Langevin step size and the number of inference steps.

We additionally conduct a wall-clock time analysis as a function of the number of particles \(K\). Table 4 shows that VarEOT training time scales with \(K\): larger \(K\) leads to longer training, but also to better FID. Importantly, VarEOT achieves superior FID compared to both LightSB and EgNOT across all considered values of \(K\), while requiring significantly less training time than EgNOT (e.g., \(131\)s vs.\(1692\)s for better quality). We attribute this to our simulation-free training objective, which avoids costly trajectory simulation required by EgNOT.

6 Discussion↩︎

Potential impact. By introducing an exact variational reformulation of the weak dual EOT objective, VarEOT makes a step toward more efficient entropic transport algorithms that avoid the key limitations of existing methods, such as simulation-based training, adversarial objectives, and restrictive parametric assumptions.

Limitations. While VarEOT enables simulation-free training, several limitations remain. First, sampling at inference time via Langevin dynamics is still required, making generation quality dependent on step size and step count. Second, the exponential terms in the objective may cause numerical instability during training. Third, the quality of optimization is inherently tied to the expressiveness of the auxiliary network \(\xi\), and the choice of its architecture remains a practical question; we note that in our experiments a standard MLP for \(\xi\) proved sufficient. Finally, the current framework is specialized to the quadratic cost (Gaussian kernel); extending to general costs is a natural direction for future work.

Acknowledgements↩︎

The work was supported by the grant for research centers in the field of AI provided by the Ministry of Economic Development of the Russian Federation in accordance with the agreement 000000C313925P4F0002 and the agreement №139-10-2025-033.

Impact Statement↩︎

This paper presents work whose goal is to advance the field of Machine Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here.

[36]

7 Proofs↩︎

Throughout the appendix, we will use the shorthand \(Z_f := Z(f,\cdot)\).

The following symmetrization bound is standard; see, e.g., [37].

We will use the vector contraction inequality of [38].

7.1 Proof of Proposition 1↩︎

7.2 Proof of Theorem 2↩︎

7.3 Proof of Theorem 3↩︎

For the reader’s convenience, we briefly recall the empirical loss \[\widehat{\mathcal{L}}(f, \xi) = \varepsilon - \frac{\varepsilon D}{2}\log(2 \pi \varepsilon)+ \frac{1}{N_1} \sum_{j=1}^{N_1}{f(x_j^1)} - \frac{\varepsilon}{N_0} \sum_{i=1}^{N_0} \xi(x_i^0) - \frac{\varepsilon}{N_0K} \sum_{i,k=1}^{N_0,K} \exp\!\left( \frac{f(x_i^0 + \sqrt{\varepsilon}\, z_{i,k})}{\varepsilon} - \xi(x_i^0) \right),\] where the samples are taken from the corresponding distributions \(x^1_j \sim p_1,\; x^0_i \sim p_0,\; z_{i,k}\sim \mathcal{N}(0,I)\).

7.4 Proof of Proposition 4↩︎

7.5 Proof of Theorem 5↩︎

7.6 Proof of Theorem 6↩︎

8 Details of the Experiments↩︎

8.1 VarEOT: Optimization and Architecture↩︎

8.1.0.1 Optimization.

Training is performed using the AdamW [41] optimizer with a learning rate \(\text{lr} = 10^{-4}\), momentum parameters \(\beta_1 = 0.7\) and \(\beta_2 = 0.8\), and weight decay set to \(10^{-4}\). We additionally employ an exponential moving average (EMA) of model parameters with momentum \(\text{ema\_momentum} = 0.999\), which is used for evaluation. All models are trained for \(10^4\) gradient steps.

8.1.0.2 Network architecture.

The transport potential is parameterized by a multilayer perceptron (MLP). The network consists of four fully connected layers with SiLU activations. The input dimensionality is denoted by \(d_{\text{in}}\) (512 for latent-space experiments), the hidden layer width is fixed to 256 units, and the output dimension is \(d_{\text{out}} = 1\).

8.1.0.3 Training setup.

All models are trained with a batch size of 256. The number of Monte Carlo samples used to estimate expectations is fixed to \(K = 256\) in the image experiments, for MSCI we use \(K = 128\).

8.1.0.4 Langevin inference.

The number of Langevin steps \(S\) (denoted as NFE) and the corresponding step size are chosen depending on the experimental setting. For all toy experiments on the SwissRoll dataset, we use \(1000\) Langevin steps with a fixed step size of \(10^{-3}\). For image-based experiments in the ALAE latent space, we consider multiple inference budgets. The step sizes for each NFE configuration are selected based on FID performance, as illustrated in Figure 9. Specifically, we use a step size of \(0.5\) for \(\text{NFE}=2\), a step size of \(0.1\) for \(\text{NFE}=10\), and a step size of \(10^{-3}\) for \(\text{NFE}=1000\).

8.1.0.5 Update schedule and practical considerations.

We perform simultaneous optimization of both the dual potential network \(\hat{f}_\theta\) and the auxiliary network \(\hat{\xi}_\psi\) at every training step, which we found to lead to stable training and strong empirical performance. To further stabilize training, we clip the difference values appearing in the exponential term of 10 from above by \(30\), and apply standard gradient clipping with norm threshold \(1\).

8.2 Baseline Methods↩︎

For completeness, we summarize the main technical details of the baseline methods used for comparison in the unpaired image-to-image translation experiments.

8.2.0.1 LightSB.

For LightSB, we follow the standard configuration reported in prior work. The model is trained using \(K = 10\) Gaussian components, a learning rate of \(\text{lr} = 10^{-3}\), and a batch size of 128. Optimization is performed for \(10^4\) gradient steps.

8.2.0.2 EgNOT.

For EgNOT, we follow the optimization/hyperparameter setup from the original work, but adopt the same neural network architecture for potential \(f\) as in section 8.1 for fair comparison with VarEOT. Training/inference hyperparameters were selected by grid search over a reasonable parameter space based on FID performance. For all ALAE experiments, Adam optimizer is used with \(\text{lr} = 5\cdot10^{-5}\), \((\beta_1, \beta_2) = (0, 0.999)\). Optimization is performed for \(10^4\) gradient steps, batch size is \(128\).

8.3 Image Data and Preprocessing↩︎

For image-based experiments, we rely on the official implementation of ALAE and the corresponding pretrained models, available at

To obtain semantic attributes for the FFHQ dataset, we use the neural network–based annotations released at

For baseline comparisons, we use the official implementations of LightSB and EgNOT, available at

and

respectively.

9 Additional Experimental Results↩︎

9.1 Wall-Clock Training Time Comparison↩︎

We report wall-clock training times for the M\(\to\)F experiment (\(\varepsilon{=}1.0\)) with \(10000\) iterations, using the same hyperparameters as in Appendix 8.1. For VarEOT, we vary the number of particles \(K \in \{32, 64, 128, 256, 512\}\) to study the speed–quality trade-off. All experiments for VarEOT and EgNOT were conducted on a single GPU, while LightSB experiments were run on CPU, as in the original paper.

9.2 Other Unpaired Image-to-Image Translation Setups↩︎

In the main text, we present qualitative results for the unpaired image-to-image translation task using the ALAE encoder in the Male\(\to\)Female (M\(\to\)F) setup. In this appendix, we provide additional qualitative examples for the same task under alternative and commonly used translation directions.

Specifically, we report results for the Woman\(\to\)Man (F\(\to\)M), Adult\(\to\)Child (A\(\to\)C), and Child\(\to\)Adult (C\(\to\)A) setups. All experiments are conducted under the same training protocol and model configuration as in the main text, and differ only in the choice of the source and target domains. For consistency, we fix the regularization parameter to \(\varepsilon = 1.0\) for all setups shown below.

Figures 5–7 demonstrate that the proposed method successfully learns meaningful transport maps across different translation directions, further confirming that the qualitative behavior observed in the M\(\to\)F setup generalizes to other unpaired image-to-image translation scenarios.

9.3 Dependence on the parameter \(\varepsilon\).↩︎

In Figure 8, we show how the solution learned by VarEOT depends on the parameter \(\varepsilon\) in the Male\(\rightarrow\)Female experiment. As expected, the diversity increases with the increase of \(\varepsilon\).

9.4 Effect of Langevin Inference Parameters↩︎

To further analyze the behavior of VarEOT at inference time, we study the sensitivity of the translation quality to the Langevin sampling parameters in the Male\(\rightarrow\)Female (M\(\to\)F) image-to-image translation task. Figure 9 reports heatmaps of FID scores as a function of the Langevin step size and the number of sampling steps, for different values of the entropic regularization parameter \(\varepsilon\).

References↩︎