Paper deep dive
Posterior sampling by source-space MCMC via prior-based few-step transport maps
Hoang Phuc Hau Luu, Marcelo Hartmann, Zhongjian Wang
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 92%
Last extracted: 10/3/2026, 2:01:27 AM
Summary
This paper introduces a framework for source-space generalized Bayesian inference that utilizes improved MeanFlow (iMF) transport maps to represent implicit priors. By performing Markov Chain Monte Carlo (MCMC) sampling in the Gaussian source space of the iMF map, the method achieves computational efficiency through few-step generation while providing theoretical guarantees on posterior stability via Wasserstein error bounds. The framework employs parallel tempering with preconditioned Crank-Nicolson updates and a hybrid split Hamiltonian Monte Carlo variant, demonstrating accuracy on synthetic data and text-guided image generation using CLIP.
Entities (7)
Relation Signals (6)
Improved MeanFlow → usedfor → Generalized Bayesian Inference
confidence 95% · We represent the prior using a one- or few-step improved MeanFlow (iMF) map and perform posterior sampling in its Gaussian source space.
CLIP → usedin → ImageNet experiments
confidence 93% · CLIP-guided ImageNet experiments demonstrate its ability to steer a pretrained iMF image prior toward text-specified preferences.
Improved MeanFlow → enables → Few-step generation
confidence 92% · iMF enables few-step generation and ensures that L differs from the marginal risk R only by an additive constant independent of theta.
Hybrid Variant → incorporates → Split Hamiltonian Monte Carlo
confidence 90% · ...introduce a hybrid variant that incorporates split Hamiltonian Monte Carlo to improve sampling efficiency.
Parallel Tempering → uses → Preconditioned Crank-Nicolson
confidence 90% · In the iMF source space, we adopt parallel tempering with preconditioned Crank–Nicolson updates...
Wasserstein Error Bounds → quantifies → Posterior Stability
confidence 88% · We establish Wasserstein error bounds between the exact and learned posteriors... We refer to this as posterior stability.
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Bayesian inference increasingly uses informative but implicit priors represented only by samples, such as historical ensembles, simulator outputs, and pretrained generative models. The same computational problem appears in the test-time guidance task (generalized Bayes), where an explicit positive weight, e.g., an exponentiated reward, tilts an implicit prior. We develop a framework for source-space generalized Bayesian inference that combines inexpensive few-step prior transports with posterior stability guarantees. Specifically, we represent the prior using a one- or few-step improved MeanFlow (iMF) map and perform posterior sampling in its Gaussian source space. We establish Wasserstein error bounds between the exact and learned posteriors in terms of the joint population iMF and auxiliary-velocity loss, decomposed into training suboptimality and model-class approximation error. In the iMF source space, we adopt parallel tempering with preconditioned Crank-Nicolson updates and introduce a hybrid variant that incorporates split Hamiltonian Monte Carlo to improve sampling efficiency. Synthetic experiments show that the proposed framework can approximate posterior distributions accurately and efficiently, while CLIP-guided ImageNet experiments demonstrate its ability to steer a pretrained iMF image prior toward text-specified preferences.
Tags
Links
- Source: https://arxiv.org/abs/2610.01034v1
- Canonical: https://arxiv.org/abs/2610.01034v1
Trouble viewing inline? Open PDF directly →
Full Text
94,484 characters extracted from source content.
Expand or collapse full text
Posterior sampling by source-space MCMC via prior-based few-step transport maps Hoang Phuc Hau Luu †thanks: Division of Mathematical Sciences, School of Physical and Mathematical Sciences, Nanyang Technological University, Singapore. hoangph.luu@ntu.edu.sg Marcelo Hartmann †thanks: Department of Computer Science, Faculty of Science, University of Helsinki, Finland. marcelo.hartmann@helsinki.fi Zhongjian Wang †thanks: Division of Mathematical Sciences, School of Physical and Mathematical Sciences, Nanyang Technological University, Singapore. zhongjian.wang@ntu.edu.sg (Corresponding) October 1, 2026 Abstract Bayesian inference increasingly uses informative but implicit priors represented only by samples, such as historical ensembles, simulator outputs, and pretrained generative models. The same computational problem appears in the test-time guidance task (generalized Bayes), where an explicit positive weight, e.g., an exponentiated reward, tilts an implicit prior. We develop a framework for source-space generalized Bayesian inference that combines inexpensive few-step prior transports with posterior stability guarantees. Specifically, we represent the prior using a one- or few-step improved MeanFlow (iMF) map and perform posterior sampling in its Gaussian source space. We establish Wasserstein error bounds between the exact and learned posteriors in terms of the joint population iMF and auxiliary-velocity loss, decomposed into training suboptimality and model-class approximation error. In the iMF source space, we adopt parallel tempering with preconditioned Crank–Nicolson updates and introduce a hybrid variant that incorporates split Hamiltonian Monte Carlo to improve sampling efficiency. Synthetic experiments show that the proposed framework can approximate posterior distributions accurately and efficiently, while CLIP-guided ImageNet experiments demonstrate its ability to steer a pretrained iMF image prior toward text-specified preferences. 1 Introduction Bayesian inference combines prior knowledge with observations to characterize uncertainty about unknown quantities (Stuart, 2010). In many applications, however, prior information is available through samples, e.g., historical ensembles, simulator outputs, or pretrained generative models rather than an explicit density (Patel and Oberai, 2019; Chung et al., 2023; Lunz et al., 2018). Such prior samples can represent complex structures that are difficult to encode analytically, but their densities and scores are not directly available for posterior Markov chain Monte-Carlo (MCMC). Formally, let x∈ℝdxx ^d_x be an unknown with prior π0 _0 available only through i.i.d. samples, and let y denote observations or a task specification. We encompass Bayesian and generalized Bayesian inference through the generalized posterior πy(x)=Zy−1ℓy(x)π0(x),Zy=∫ℓy(x)π0(x),π^y(dx)=Z_y^-1 _y(x) _0(dx), Z_y= _y(x) _0(dx), (1) where ℓy _y is an analytically given nonnegative weight and 0<Zy<∞0<Z_y<∞. Standard Bayesian inference takes ℓy(x)=p(y|x) _y(x)=p(y|x), while generalized Bayes (Bissiri et al., 2016) uses ℓy(x)=exp(−τLy(x)) _y(x)= (-τ L_y(x)) for a task loss LyL_y and τ>0τ>0. Reward guidance (Korbak et al., 2022) fits the same formulation with ℓy(x)=exp(τRy(x)) _y(x)= (τ R_y(x)), favoring samples with larger reward RyR_y. For example, y may be a text prompt, with Ry(x)=CLIP(x,y)R_y(x)=CLIP(x,y) measuring the cosine similarity between the CLIP embeddings of image x and prompt y (Radford et al., 2021). Henceforth, we use ℓy _y for either a likelihood or a general weight and refer to πyπ^y simply as the posterior. An approach is to learn the data-space prior score by denoising score matching (DSM) (Vincent, 2011) and combine it with the gradient of an analytic log-likelihood or reward to construct Langevin or Hamiltonian dynamics. These samplers are typically unadjusted because a learned score alone does not provide the density ratios needed for Metropolis correction. Moreover, sharp variation in the score near complex data manifolds can produce stiff dynamics, requiring small integration steps. Another approach is to first learn a deterministic generator Sθ:ℝdz→ℝdxS_θ:R^d_z ^d_x satisfying (Sθ)γ#≈π0(S_θ)_\#γ≈ _0, where γ=(0,I)γ=N(0,I) denotes the source distribution, ℝdzR^d_z is called the source space and # denotes the pushforward. Under the learned model SθS_θ, the complex data-space prior is represented by the simple Gaussian reference measure γ. The corresponding learned posterior in the source space is νθy(dz)∝ℓy(Sθ(z))γ(dz). _θ^y(dz) _y(S_θ(z))γ(dz). (2) We sample from νθy _θ^y using MCMC and map the resulting samples to data space through SθS_θ. In contrast to the conventional data-space approach, νθy _θ^y is available up to normalization, enabling Metropolis correction and thus an MCMC kernel that leaves νθy _θ^y invariant. Furthermore, the source-space posterior is expected to have more favorable geometry for MCMC: its prior is Gaussian, while SθS_θ has already absorbed much of the curvature and anisotropy of π0 _0. Nevertheless, the likelihood/weight ℓy∘Sθ _y S_θ may reintroduce curvature, stiffness, anisotropy, or multimodality, particularly when the observation model is highly informative or nonlinear. Early work by Patel et al. (2022) used a Wasserstein GAN to construct SθS_θ and applied Hamiltonian Monte Carlo (HMC) with the No-U-Turn Sampler (Hoffman et al., 2014) to the resulting low-dimensional latent posterior. However, learning SθS_θ with a WGAN requires solving an adversarial min-max optimization problem, which can be challenging to stabilize. Flow-based models (Chen et al., 2018) offer an attractive alternative: flow matching learns a time-dependent vector field through a standard regression objective (Lipman et al., 2023; Tong et al., 2024). This provides a scalable approach to learning expressive transport maps for complex and multimodal priors. Under suitable regularity conditions on the vector field, the resulting ODE flow also defines a diffeomorphism between source and data spaces (Papamakarios et al., 2021). Source-space guidance for flow-matching models has been explored through both optimization and sampling. D-Flow (Ben-Hamu et al., 2024) optimizes the source variable, while Source-Guided Flow Matching (SGFM) (Wang et al., 2026b) samples from the source posterior distribution (2). Source Parallel Tempering (SPT) (Wang et al., 2026a) combines local preconditioned Crank-Nicolson (pCN) (Beskos et al., 2008) updates with replica exchanges (Earl and Deem, 2005) across tempered source distributions to address multimodality. These source-space methods are instances of test-time guidance, incorporating user-specified likelihoods or weights during inference. For methods that represent the prior through a learned transport (Wang et al., 2026b; Wang et al., 2026a), a central theoretical question is how prior-learning error controls the discrepancy between the learned posterior (Sθ)#νθy(S_θ)_\# _θ^y and the true posterior πyπ^y in (1). We refer to this as posterior stability. Existing works have tried to address this question, but a gap remains. Wang et al. (2026b) assume a uniform vector-field approximation bound, ‖vθ(x,t)−v(x,t)‖∞≤ϵ\|v_θ(x,t)-v(x,t)\|_∞≤ε where v(x,t)v(x,t) is the true velocity field and vθ(x,t)v_θ(x,t) is the learned velocity field. However, the flow-matching objective controls an average squared error rather than a uniform error. Wang et al. (2026a) address this mismatch by establishing posterior stability in terms of the L2L^2 transport-map error, i.e., z∼γ‖S(z)−Sθ(z)‖2≤ϵ2E_z γ\|S(z)-S_θ(z)\|^2≤ε^2, where S is the exact transport map induced by v(x,t)v(x,t). Their connection to the flow matching objective takes SθS_θ to be the map induced by the exactly integrated learned flow, via Grönwall’s inequality. In practice, however, numerical ODE solvers are employed, introducing an additional discretization error. Under standard regularity and stability assumptions, a numerical ODE solver of order k has global error O(hk)O(h^k) (Hairer et al., 1993). The training objective alone does not control this error. Reducing h requires more integration steps, increasing the cost of source-space MCMC. Minibatch optimal-transport coupling (Tong et al., 2024), as used by Wang et al. (2026b), can encourage straighter trajectories and facilitate few-step integration, but does not by itself provide a bound on this numerical error; their theoretical bound likewise concerns the exactly integrated learned flow. Thus, guarantees for the exactly integrated learned flow do not by themselves establish posterior stability for a fixed few-step implementation. Likewise, Purohit et al. (2025) use few-step consistency models but do not establish how the consistency-distillation training objective controls the resulting posterior error. In this work, we pursue two complementary objectives: making source-space inference computationally efficient through few-step transports and establishing posterior stability guarantees for the implemented maps via the true population training objective. To this end, we use improved MeanFlow (iMF) (Geng et al., 2026b), which learns interval-average velocities without distillation and natively enables one-step or few-step generation without numerical ODE integration. We further explain why the original MeanFlow objective (Geng et al., 2026a) does not directly yield the same guarantee, highlighting that native few-step generation alone can fall short of establishing such stability results. Practically, the iMF maps make source-space gradients inexpensive to evaluate when available, allowing MCMC kernels to optionally exploit them to improve mixing. Contributions We use iMF transport mappings to represent complex priors available only through samples in generalized Bayesian inference settings. Under stated regularity assumptions, we connect population training suboptimality and model-class approximation error for the joint iMF and auxiliary-velocity objective to posterior stability for the implemented one-step or few-step maps. Building on these transports, we adopt SPT with pCN updates (Wang et al., 2026a) to facilitate exploration across modes while exploiting the Gaussian reference structure of the high-dimensional source space. When likelihood or reward gradients are available, we additionally introduce Gaussian-reference split HMC (Beskos et al., 2011) at the coldest temperature of SPT to exploit local posterior geometry and potentially improve mixing. Both local kernels preserve the Gaussian reference measure in their reference dynamics, making them natural choices for inference in high-dimensional source spaces. We evaluate the proposed framework on synthetic posterior sampling and CLIP-guided ImageNet generation, demonstrating accurate and efficient posterior approximation in the synthetic settings and the ability to steer a pretrained image prior toward text-specified preferences. Other related work Transport MCMC Transport maps have long been used to improve MCMC by transforming difficult posterior target distributions into coordinates easier for sampling (Parno and Marzouk, 2018; Zhang et al., 2025). NeuTra HMC learns a variational transport tailored to posterior geometry, with an amortized extension that accommodates different observations within a specified probabilistic model (Hoffman et al., 2019). Score-based and diffusion inverse solvers Diffusion-based inverse solvers offer another way to combine a learned sample-only prior with observations, typically by incorporating likelihood information into the reverse generative process. Diffusion posterior sampling (DPS) (Chung et al., 2023) is a prominent example, enabling posterior-oriented sampling. Its approximate guidance, however, does not generally preserve the desired posterior and can introduce systematic bias beyond that arising from prior learning and numerical discretization. Denoising-score MCMC (Ramzi et al., 2020) learns the prior score from samples using denoising score matching (Vincent, 2011) and combines it with an analytically available likelihood score to perform MCMC directly in data space. Training-time posterior learning. Complementary to test-time guidance, amortized posterior learning trains observation-conditional generators to directly approximate posterior distributions. Neural posterior estimation and related methods learn from simulated parameter-observation pairs, enabling rapid inference for new observations (Greenberg et al., 2019; Radev et al., 2023; Cheng et al., 2026). The training cost is shared across subsequent inference tasks, but changes to the observation model beyond those covered during training may require retraining or adaptation. 2 Background 2.1 Flow matching Flow matching (FM) (Lipman et al., 2023) learns a velocity field approximating the marginal ODE w˙t=v(wt,t) w_t=v(w_t,t), which transports a simple source distribution psrcp_src, typically Gaussian, to a target distribution ptgtp_tgt available through samples. Let X∼ptgtX p_tgt and Z∼psrcZ p_src be independent target and source samples, and define Wt=(1−t)X+tZW_t=(1-t)X+tZ and the samplewise interpolation velocity V=Z−XV=Z-X. FM approximates the marginal instantaneous velocity v(w,t)=[V∣Wt=w]v(w,t)=E[V W_t=w] with a learned field vθv_θ, whose approximation error is measured by the marginal FM risk ℛFM(θ)=‖vθ(WT,T)−v(WT,T)‖22. _FM(θ)=E\|v_θ(W_T,T)-v(W_T,T)\|^2_2. (3) where T∼Unif[0,1]T [0,1] independent of X,ZX,Z and the above expectation is taken over T,X,ZT,X,Z. Since v is intractable, we replace v by V, resulting in the population training loss ℒFM(θ)=‖vθ(WT,T)−V‖22.L_FM(θ)=E \|v_θ(W_T,T)-V \|_2^2. (4) Minimizing ℒFML_FM is equivalent to minimizing ℛFMR_FM since ℒFM(θ)=ℛFM(θ)+CL_FM(θ)=R_FM(θ)+C where C=‖V−v(WT,T)‖22C=E\|V-v(W_T,T)\|_2^2 does not depend on θ (Lipman et al., 2023). Once trained, generation integrates w˙t=vθ(wt,t) w_t=v_θ(w_t,t) from t=1t=1 to t=0t=0, using a numerical ODE solver. 2.2 MeanFlow and improved MeanFlow MeanFlow (Geng et al., 2026a) instead parameterizes the average velocity u(wt,r,t)=1t−r∫rtv(wτ,τ)τ.u(w_t,r,t)= 1t-r _r^tv(w_τ,τ)\,dτ. Fundamental calculus yields the so-called MeanFlow identity u(wt,r,t)=v(wt,t)−(t−r)(∂tu+(∂wu)v(wt,t)). u(w_t,r,t)=v(w_t,t)-(t-r) ( _tu+( _wu)v(w_t,t) ). (5) Based on this identity, the original MeanFlow objective regresses an average velocity uθu_θ to the right hand side of (5), and since ∂tu _tu and ∂wu _wu are unknown, it replaces those by ∂tuθ _tu_θ and ∂wuθ _wu_θ to get ℛMF(θ)=‖uθ(WT,R,T)−[v(WT,T)−(T−R)(∂tuθ+(∂wuθ)v(WT,T))]‖22,R_MF(θ)=E \|u_θ(W_T,R,T)-\! [v(W_T,T)-(T-R) ( _tu_θ+( _wu_θ)v(W_T,T) ) ] \|_2^2, (6) where the expectation is taken over R,T,Z,XR,T,Z,X (recall Wt=(1−t)X+tZW_t=(1-t)X+tZ). We omit the stop-gradient operator when defining population risks, since they leave objective values unchanged. Replacing v by V gives a tractable training loss, ℒMF(θ)=‖uθ(WT,R,T)−[V−(T−R)(∂tuθ+(∂wuθ)V)]‖22.L_MF(θ)=E \|u_θ(W_T,R,T)- [V-(T-R) ( _tu_θ+( _wu_θ)V ) ] \|_2^2. (7) Unlike flow matching, the gap ℒMF(θ)−ℛMF(θ)L_MF(θ)-R_MF(θ) generally depends on θ, making them inequivalent as scalar objectives11 1 The stop-gradient makes the two objectives share expected gradients (Lu and Wang, 2026), but this is an optimization device; our theory concerns objective values, not training dynamics.. Moreover, the training objective (7) is not standard since its effective target depends on the network itself. To address this, improved MeanFlow (iMF) (Geng et al., 2026b) replaces the interpolation velocity in (∂wuθ)V( _wu_θ)V in ℒMFL_MF by a learnable marginal velocity vθv_θ, defining ℒiMF(θ) _iMF(θ) =‖uθ(WT,R,T)+(T−R)[∂tuθ+(∂wuθ)vθ(WT,T)]−V‖22. =E \|u_θ(W_T,R,T)+(T-R)\! [ _tu_θ+( _wu_θ)v_θ(W_T,T) ]-V \|_2^2. (8) The corresponding marginal iMF risk is then ℛiMF(θ)=‖uθ(WT,R,T)+(T−R)[∂tuθ+(∂wuθ)vθ(WT,T)]−v(WT,T)‖22 _iMF(θ)=E \|u_θ(W_T,R,T)+(T-R)\! [ _tu_θ+( _wu_θ)v_θ(W_T,T) ]-v(W_T,T) \|_2^2 (9) Importantly, ℒiMF(θ)−ℛiMF(θ)L_iMF(θ)-R_iMF(θ) is independent of θ (Appendix A.1), recovering the equivalence between these objectives. vθ(wt,t)v_θ(w_t,t) can be represented by the boundary value uθ(wt,t,t)u_θ(w_t,t,t) or an auxiliary velocity head. In this work, we consider the latter, and the auxiliary head is explicitly supervised through FM loss (4), ℒFM(θ)=‖vθ(WT,T)−V‖22L_FM(θ)=E \|v_θ(W_T,T)-V \|_2^2, yielding the joint training loss ℒ(θ)=ℒiMF(θ)+λℒFM(θ),λ>0. (θ)=L_iMF(θ)+ _FM(θ), λ>0. (10) We also denote the joint risk ℛ(θ)=ℛiMF(θ)+λℛFM(θ). (θ)=R_iMF(θ)+ _FM(θ). (11) The auxiliary velocity field vθv_θ is used only during training; sampling relies solely on uθu_θ, with the one-step map given by Sθ(z)=z−uθ(z,0,1),z∼psrc,S_θ(z)=z-u_θ(z,0,1), z p_src, (12) More generally, a fixed schedule 1=tK>tK−1>⋯>t0=01=t_K>t_K-1>·s>t_0=0 defines a few-step map through wK=z,wk−1=wk−(tk−tk−1)uθ(wk,tk−1,tk),k=K,…,1,Sθ(z)=w0.w_K=z, w_k-1=w_k-(t_k-t_k-1)u_θ(w_k,t_k-1,t_k), k=K,…,1, S_θ(z)=w_0. (13) 3 Source-space sampling via improved MeanFlow 3.1 Bayesian inference through a learned prior transport Let π0 _0 be a prior accessed through samples and ℓy _y a nonnegative likelihood or weight associated with observations or a task y. We aim to sample from the posterior πy(dx)∝ℓy(x)π0(dx)π^y(dx) _y(x) _0(dx). The idea of source-space sampling (Purohit et al., 2025; Wang et al., 2026b) is to learn a transport map (Sθ)#γ≈π0(S_θ)_\#γ≈ _0 where γ=(0,I)γ=N(0,I) and perform posterior sampling in the source space of γ, targeting the pull-back posterior νθy(dz)∝ℓy(Sθ(z))γ(dz) _θ^y(dz) _y(S_θ(z))γ(dz) then transport these source samples to the data space via SθS_θ. In case the transport map is exact, i.e., (Sθ)#γ=π0(S_θ)_\#γ= _0, an exact posterior source sample z∼νθy(dz)z _θ^y(dz) will give an exact sample, Sθ(z)S_θ(z), from πy(dx)π^y(dx) (Wang et al., 2026b, Thm. 1). It is then critical to learn a good prior transport SθS_θ and perform efficient sampling in the source space of SθS_θ. When SθS_θ is not perfect, we seek to quantify the resulting posterior error. Let S denote the exact transport, i.e., S#γ=π0S_\#γ= _0, associated with the generative construction, assumed to exist under suitable regularity conditions. Define the ideal source posterior by νy(dz)∝ℓy(S(z))γ(dz)ν^y(dz) _y(S(z))γ(dz). Our goal is to control the source-space discrepancy W2(νθy,νy)W_2( _θ^y,ν^y) and the data-space discrepancy W2((Sθ)#νθy,πy)W_2((S_θ)_\# _θ^y,π^y) in terms of the training error associated with learning SθS_θ (for iMF, it links to ℛ(θ)R(θ), see next section). Here W2W_2 denotes the Wasserstein distance (Villani and others, 2009). 3.2 Improved MeanFlow source and its posterior stability We parameterize SθS_θ using improved MeanFlow (12) and (13), trained with the joint training loss ℒ(θ)=ℒiMF(θ)+λℒFM(θ)L(θ)=L_iMF(θ)+ _FM(θ) in (10). iMF enables few-step generation and ensures that ℒL differs from the marginal risk ℛR (11) only by an additive constant independent of θ. Consequently, for a fixed model class Θ , ℛ(θ)=ℒ(θ)−infϑ∈Θℒ(ϑ)⏟:=ϵopt,population training suboptimality+infϑ∈Θℛ(ϑ)⏟:=ϵapp,model-class approximation error:=δ. (θ)= L(θ)- _ ∈ L( )_:= _opt,~population training suboptimality+ _ ∈ R( )_:= _app,~model-class approximation error:=δ. (14) See App. A.2 for proof. The decomposition (14) shows that the marginal risk ℛ(θ)R(θ) is controlled by population suboptimality of the training objective ℒL used in practice and model-class approximation error. By the MeanFlow identity (5), the exact instantaneous and interval-average velocity fields jointly achieve zero marginal risk (Geng et al., 2026b), so a sufficiently expressive model class can in principle make ϵapp _app small. We therefore seek to bound posterior error in terms of ℛ(θ)R(θ). We also note that, for original MeanFlow (Geng et al., 2026a), the parameter-dependent gap ℒMF(θ)−ℛMF(θ)L_MF(θ)-R_MF(θ) prevents the decomposition (14) in general. Small training suboptimality and approximation error therefore need not imply small marginal risk. We establish posterior stability by first controlling the error of the implemented iMF transport, building on the flow-map error analysis of (Boffi et al., 2024). For the analysis, we use independent coupling X∼π0,Z∼γX _0,Z γ and time R,T∼iidUnif[0,1]R,T iid Unif[0,1]22 2 A more general schedule is possible, but we use a simpler one for clarity of presentation. . Also, let V=Z−XV=Z-X, Wt=(1−t)X+tZW_t=(1-t)X+tZ, v(w,t)=[V|Wt=w]v(w,t)=E[V|W_t=w], and pt=Law(Wt)p_t=Law(W_t). Thus, p0=π0p_0= _0 and p1=γp_1=γ. Assumption 3.1 (Prior and marginal flow). The prior π0∈2(ℝd) _0 _2(R^d). For each t∈[0,1]t∈[0,1] and ptp_t-a.e. w, the ODE dsξst,w=v(ξst,w,s),ξt,w=w dds _s^t,w=v( _s^t,w,s),~ _t^t,w=w, admits absolutely continuous solutions on [0,t][0,t] whose solution maps φs,t(w):=ξst,w _s,t(w):= _s^t,w form a measurable flow satisfying (φs,t)#pt=ps,0≤s≤t≤1( _s,t)_\#p_t=p_s,0≤ s≤ t≤ 1. Under Assumption 3.1, the map S:=φ0,1S:= _0,1 is an exact prior transport: S#γ=π0S_\#γ= _0. Assumption 3.2 (Regularity of the learned average velocity). The field uθu_θ is continuously differentiable on ℝd×(0,1)2R^d×(0,1)^2, and both uθu_θ and its first derivatives extend continuously to ℝd×[0,1]2R^d×[0,1]^2. For r,t∈[0,1]r,t∈[0,1], define the marginal iMF regression residual q¯θ(w,r,t) q_θ(w,r,t) :=uθ(w,r,t)−v(w,t)+(t−r)[∂tuθ(w,r,t)+∂wuθ(w,r,t)vθ(w,t)] :=u_θ(w,r,t)-v(w,t)+(t-r) [ _tu_θ(w,r,t)+ _wu_θ(w,r,t)\,v_θ(w,t) ] Gθ(r) G_θ(r) :=∫r1‖q¯θ(⋅,r,t)‖L2(pt)2t. := _r^1\| q_θ(·,r,t)\|_L^2(p_t)^2\,dt. Thus, the iMF risk in (9) is rewritten as ℛiMF(θ)=‖q¯θ(WT,R,T)‖2R_iMF(θ)=E\| q_θ(W_T,R,T)\|^2. We impose the following. Assumption 3.3 (Regularity across starting times). Gθ:[0,1]→[0,∞)G_θ:[0,1]→[0,∞) is LGL_G-Lipschitz, i.e., |Gθ(r)−Gθ(r′)|≤LG|r−r′|,∀r,r′∈[0,1].|G_θ(r)-G_θ(r )|≤ L_G|r-r |, ∀ r,r ∈[0,1]. Remark 3.4. Assumption 3.3 holds if q¯θ q_θ is continuously differentiable in r and both q¯θ q_θ and ∂rq¯θ _r q_θ are uniformly bounded in L2(pt)L^2(p_t) over 0≤r≤t≤10≤ r≤ t≤ 1. Under the finite-second-moment assumption above, these conditions hold for finite affine-SiLU networks with direct spatial and time inputs; the resulting constants may depend on the network parameters. See App. A.3 for details and proof. Fix a generation schedule 1=tK>tK−1>⋯>t0=01=t_K>t_K-1>·s>t_0=0, put hk=tk−tk−1h_k=t_k-t_k-1. Assumption 3.5. ∑k=1K∫tk−1tk(t−tk−1)2‖(∂wuθ)(⋅,tk−1,t)(v−vθ)(⋅,t)‖L2(pt)2t≤Bu2ℛFM(θ). _k=1^K _t_k-1^t_k(t-t_k-1)^2 \|( _wu_θ)(·,t_k-1,t)(v-v_θ)(·,t) \|_L^2(p_t)^2\,dt≤ B_u^2R_FM(θ). (15) where ℛFMR_FM is the marginal flow-matching risk (3) of the auxiliary velocity head. Remark 3.6. A sufficient condition for (15) is supw, 0≤r<t≤1(t−r)‖∂wuθ(w,r,t)‖op≤Bu _w,\,0≤ r<t≤ 1(t-r)\| _wu_θ(w,r,t)\|_op≤ B_u. In particular, if uθu_θ is LuL_u-Lipschitz in w, uniformly over (r,t)(r,t) (e.g., affine–SiLU networks), then this condition holds with Bu=LuB_u=L_u, since t−r≤1t-r≤ 1. The weighted bound in Assumption 3.5 is weaker: it only controls the Jacobian acting on v−vθv-v_θ in the L2(pt)L^2(p_t) sense. See Appendix A.4 for more details. Theorem 3.7. Let Fk,θ(w)=w−hkuθ(w,tk−1,tk)F_k,θ(w)=w-h_ku_θ(w,t_k-1,t_k) be LkL_k-Lipschitz, and let Sθ=F1,θ∘⋯∘FK,θS_θ=F_1,θ ·s F_K,θ be the few-step iMF map (13). Define Λk:=∏j<kLj,AK:=∑k=1KΛkhk,BK:=(∑k=1KΛk2hk)1/2,Λ1=1 _k:= _j<kL_j,A_K:= _k=1^K _k h_k,B_K:= ( _k=1^K _k^2h_k )^1/2, _1=1. Under Assumptions 3.1, 3.2, 3.3 and 3.5, if δ≤LG/2δ≤ L_G/2, γ‖Sθ(Z)−S(Z)‖2≤[AK(2LGδ)1/4+BuBKλδ1/2]2.E_γ\|S_θ(Z)-S(Z)\|^2≤ [A_K(2L_Gδ)^1/4+ B_uB_K λδ^1/2 ]^2. (16) If uθu_θ is uniformly LuL_u-Lipschitz in space, then Λk≤eLu _k≤ e^L_u, AK≤eLuKA_K≤ e^L_u K and BK≤eLuB_K≤ e^L_u. Proof of Theorem 3.7 is given in Appendix A.5. The bound in (16) links transport error to population training suboptimality and model-class approximation error through δ in (14). The bound vanishes when δ=0δ=0, without requiring the step sizes hkh_k to tend to zero. This shows the few-step nature of the result. In particular, K=1K=1 gives A1=B1=1A_1=B_1=1 and a transport guarantee for the one-step map (12), with no Lipschitz condition on SθS_θ needed in (16). We support Theorem 3.7 with an experiment in Appendix B.1, showing that decreasing δ is accompanied by decreasing transport-map error γ‖Sθ(Z)−S(Z)‖2E_γ\|S_θ(Z)-S(Z)\|^2 for both one- and six-step generation. Define terms in Theorem 3.7: ℰ:=Z∼γ‖Sθ(Z)−S(Z)‖2,bK(δ):=AK(2LGδ)1/4+BuBKλδ1/2.E:=E_Z γ\|S_θ(Z)-S(Z)\|^2, b_K(δ):=A_K(2L_Gδ)^1/4+ B_uB_K λδ^1/2. (17) We now connect the posterior error to the transport error. Fix y, we make the following assumption. Assumption 3.8 (Weight regularity). There exist finite constants My>0M_y>0 and Ly≥0L_y≥ 0 such that 0≤ℓy(x)≤My,|ℓy(x)−ℓy(x′)|≤Ly∥x−x′∥,∀x,x′∈ℝd.0≤ _y(x)≤ M_y, | _y(x)- _y(x )|≤ L_y\|x-x \|, ∀ x,x ^d. Moreover, Zy:=∫ℓy(x)π0(x)>0Z_y:= _y(x)\, _0(dx)>0. With the notation of Section 3.1, the normalized source posteriors are νy(z)=ℓy(S(z))Zyγ(z),νθy(z)=ℓy(Sθ(z))Zθ,yγ(z),Zθ,y:=∫ℓy(Sθ(z))γ(z),ν^y(dz)= _y(S(z))Z_yγ(dz), _θ^y(dz)= _y(S_θ(z))Z_θ,yγ(dz), Z_θ,y:= _y(S_θ(z))\,γ(dz), where νθy _θ^y is defined whenever Zθ,y>0Z_θ,y>0. Here S#νy=πyS_\#ν^y=π^y. When ℓy>0 _y>0, write Φy=−logℓy _y=- _y. Theorem 3.9. Under Assumption 3.8, if Lyℰ1/2≤Zy/2L_yE^1/2≤ Z_y/2, then Zθ,y≥Zy/2>0Z_θ,y≥ Z_y/2>0 and W2(νθy,νy)≤Cy(d)1/2ℰ1/4,Cy(d):=4LyZy[d2+4d+1+My(1+d)Zy].W_2( _θ^y,ν^y)≤ C_y(d)^1/2E^1/4, C_y(d):= 4L_yZ_y [ d^2+4d+1+ M_y(1+d)Z_y ]. (18) If, additionally, SθS_θ is LSθL_S_θ-Lipschitz, then W2((Sθ)#νθy,πy)≤LSθCy(d)1/2ℰ1/4+MyZyℰ1/2W_2 ((S_θ)_\# _θ^y,π^y )≤ L_S_θC_y(d)^1/2E^1/4+ M_yZ_yE^1/2. Proof is in Appendix A.6. Combining these results, we derive the following. Corollary 3.10. Suppose the hypotheses of Theorem 3.7 and Assumption 3.8 hold. Let LSθL_S_θ be a Lipschitz constant of SθS_θ, and let Cy(d)C_y(d) be as in Theorem 3.9. If δ≤LG/2,LybK(δ)≤Zy/2δ≤ L_G/2,L_yb_K(δ)≤ Z_y/2, then Zθ,y≥Zy/2>0Z_θ,y≥ Z_y/2>0 and W2(νθy,νy) W_2( _θ^y,ν^y) ≤Cy(d)1/2bK(δ)1/2, ≤ C_y(d)^1/2b_K(δ)^1/2, (19) W2((Sθ)#νθy,πy) W_2 ((S_θ)_\# _θ^y,π^y ) ≤LSθCy(d)1/2bK(δ)1/2+MyZybK(δ). ≤ L_S_θC_y(d)^1/2b_K(δ)^1/2+ M_yZ_yb_K(δ). (20) 3.3 Source parallel tempering with improved MeanFlow Given the learned iMF, we aim to sample from νθy _θ^y. Although the source reference is Gaussian, the pulled-back likelihood can induce multimodality, and the source retains the data dimension. We use source parallel tempering (Wang et al., 2026a) to improve mixing between modes, together with local kernels like pCN and Gaussian-reference split HMC to address high dimensionality. Both admit formulations directly on infinite-dimensional Hilbert spaces, motivating their use for high-dimensional Gaussian-reference targets (Beskos et al., 2008; Cotter et al., 2013; Beskos et al., 2011). Parallel tempering (Geyer, 1991) runs replicas at inverse temperatures 0=β0<⋯<βJ=10= _0<·s< _J=1, targeting νθ,βy(dz)∝ℓy(Sθ(z))βγ(dz) _θ,β^y(dz) _y(S_θ(z))^βγ(dz). Local MCMC kernels explore each target, while Metropolis swaps between adjacent replicas promote transitions between modes. Each sampling sweep locally updates every replica once, then attempts one Metropolis swap per adjacent temperature pair. The coldest replica (β=1β=1) targets νθy _θ^y. The two variants below differ in their local MCMC kernels. iMF+SPT+pCN. We use pCN at every replica as in (Wang et al., 2026a). Writing Φy=−logℓy _y=- _y and Ψθ=Φy∘Sθ _θ= _y S_θ, the proposal is z′=cos(hβ)z+sin(hβ)ξz = (h_β)z+ (h_β)ξ, where ξ∼γξ γ and hβ∈(0,π/2]h_β∈(0,π/2], with acceptance probability 1∧exp−β[Ψθ(z′)−Ψθ(z)]1 \-β[ _θ(z )- _θ(z)]\. At β=0β=0, every proposal is accepted; and choosing h0=π/2h_0=π/2 gives an independent Gaussian draw z′=ξ∼γz =ξ γ. iMF+SPT+hybrid. When Ψθ _θ is differentiable, we replace the coldest pCN update (β=1β=1) with split HMC (Beskos et al., 2011) to exploit local posterior geometry, while retaining pCN at hotter replicas to limit computational cost. Refresh p∼γp γ and define the Hamiltonian H(z,p)=12‖z‖2+12‖p‖2+Ψθ(z)H(z,p)= 12\|z\|^2+ 12\|p\|^2+ _θ(z). Starting from (z0,p0)=(z,p)(z_0,p_0)=(z,p), apply L steps of Strang’s splitting: pℓ+1/2 p_ +1/2 =pℓ−h2∇Ψθ(zℓ),⊳half momentum kick =p_ - h2∇ _θ(z_ ), ~half momentum kick (21) (zℓ+1p~ℓ+1) pmatrixz_ +1\\ p_ +1 pmatrix =(cos(h)Isin(h)I−sin(h)Icos(h)I)(zℓpℓ+1/2),⊳exact harmonic rotation = pmatrix (h)I& (h)I\\ - (h)I& (h)I pmatrix pmatrixz_ \\ p_ +1/2 pmatrix, ~exact harmonic rotation pℓ+1 p_ +1 =p~ℓ+1−h2∇Ψθ(zℓ+1),⊳half momentum kick = p_ +1- h2∇ _θ(z_ +1), ~half momentum kick Accept (z′,p′)=(zL,−pL)(z ,p )=(z_L,-p_L) with probability 1∧expH(z,p)−H(z′,p′)1 \H(z,p)-H(z ,p )\; otherwise retain z. Sampler geometric convergence. Wang et al. (2026a) establish geometric W2W_2 convergence to the learned posterior for SPT with all pCN updates. We extend this guarantee to SPT+hybrid. Let νn=Law(zn) _n=Law(z_n), where znz_n is the coldest replica’s state after n sampling sweeps. Under the conditions in Appendix A.7, including a source weight bounded above and away from zero, the cold marginal of either sampler satisfies W2(νn,νθy)≤CρnW_2( _n, _θ^y)≤ Cρ^n for some C<∞C<∞ and ρ∈(0,1)ρ∈(0,1). Together with posterior stability (Theorem 3.9, Corollary 3.10), this yields an end-to-end error bound with a geometrically decaying sampling term and a transport term that vanishes with δ; see Corollary A.3. 4 Experiments We evaluate the proposed samplers on synthetic inference problems and CLIP-guided ImageNet generation. Experiments were conducted on LUMI33 3 https://lumi-supercomputer.eu/ using one Slurm-allocated GPU device per run, with AMD Instinct MI250X hardware, eight CPU cores, and 64 GB of host memory. Our code will be made available soon. 4.1 Synthetic experiments Problem setting. We consider two priors on x∈ℝ32x ^32, each formed from 16 independent two-dimensional blocks. For uj,vj∼iid(0,1)u_j,v_j iid N(0,1), the priors are x2j−1=uj,x2j=0.62(uj2−1)+0.30vj,Banana,1.05sin(1.65uj)+0.15vj,Sine,j=1,…,16.x_2j-1=u_j, x_2j= cases0.62(u_j^2-1)+0.30v_j,&Banana,\\ 1.05 (1.65u_j)+0.15v_j,&Sine, cases j=1,…,16. Both problems use linear observations with independent Gaussian noise, giving the likelihood p(y∣x)=∏j=116(yj,x2j+0.35x2j−1, 0.22).p(y x)= _j=1^16N\! (y_j;\,x_2j+0.35x_2j-1,\,0.2^2 ). We condition on the fixed observations y1=y2=2y_1=y_2=2 and yj=−1y_j=-1 for j=3,…,16j=3,…,16. Prior learning and sampling. We evaluate iMF+SPT+pCN and iMF+SPT+hybrid against two baselines: FM+SPT+pCN and diffusion posterior sampling (DPS) (Chung et al., 2023). The two iMF methods share a six-step transport. The FM baseline uses a 100-step RK4 transport. For each problem, the iMF, FM, and diffusion prior models are trained on the same fixed dataset of 20,000 prior samples. Each method returns 6,000 samples. Figure 4 compares samples from the learned priors with those from the true banana and sine priors. See Appendix B.2 for more details. Results. The blockwise structure of the problems enables accurate numerical posterior references for assessing approximation quality. Table 1 summarizes the results, and Figures 5 and 6 compare posterior sample histograms from each method with numerical references for the true banana and sine. Using the same learned iMF transport, the hybrid sampler increases mean coordinate ESS (effective sample size) by 3.24× on banana and 4.44× on sine, and coordinate ESS per second by 2.25× and 3.13×, respectively, compared with iMF+SPT+pCN. It also achieves the lowest mean sliced Wasserstein and binned TV errors, although the differences in sliced Wasserstein error are small relative to the variability across seeds. With its configured 100-step RK4 transport, FM+SPT+pCN requires 33.7× and 34.5× the sampling time of the iMF hybrid on banana and sine, respectively. DPS exhibits substantially larger distributional errors, with systematic deviations from the reference posterior visible in the histograms. Table 1: Synthetic posterior sampling. See Appendix B.2 for more details. Method SW2↓SW_2 Bin TV ↓ Coord. ESS ↑ Coord. ESS/s ↑ Time (s) ↓ Banana DPS 0.462±0.0710.462± 0.071 0.311±0.0020.311± 0.002 — — 38.4±1.038.4± 1.0 FM+SPT+pCN 0.219±0.0810.219± 0.081 0.106±0.0080.106± 0.008 372±16372± 16 0.052±0.0020.052± 0.002 7087±457087± 45 iMF+SPT+pCN 0.154±0.0120.154± 0.012 0.109±0.0060.109± 0.006 346±4346± 4 2.37±0.122.37± 0.12 146.1±6.0146.1± 6.0 iMF+SPT+hybrid 0.142±0.0310.142± 0.031 0.085±0.0050.085± 0.005 ±1122± 48 5.33±0.085.33± 0.08 210.4±5.8210.4± 5.8 Sine DPS 0.931±0.0250.931± 0.025 0.286±0.0010.286± 0.001 — — 38.9±0.238.9± 0.2 FM+SPT+pCN 0.104±0.0210.104± 0.021 0.129±0.0140.129± 0.014 395±16395± 16 0.0526±0.00230.0526± 0.0023 7505±787505± 78 iMF+SPT+pCN 0.093±0.0160.093± 0.016 0.131±0.0100.131± 0.010 378±35378± 35 2.47±0.322.47± 0.32 153.3±5.8153.3± 5.8 iMF+SPT+hybrid 0.088±0.0110.088± 0.011 0.127±0.0040.127± 0.004 ±1679± 47 7.73±0.317.73± 0.31 217.3±5.4217.3± 5.4 4.2 CLIP-guided ImageNet generation We evaluate CLIP-guided 256×256256× 256 class-conditional ImageNet generation: for a fixed ImageNet class (e.g., bald eagle) and text prompt, the goal is to sample images from that class that favor high image-text alignment. We use the cosine similarity from a frozen pretrained CLIP model (Radford et al., 2021) as the guidance reward. Let z∼(0,I)z (0,I) and let Gc(z)G_c(z) denote a frozen pretrained transport from the Gaussian source distribution to images in class c. For a text prompt p, the learned source posterior is νθp(z)∝exp(−12‖z‖2+τCLIP(Gc(z),p))ν^p_θ(z) \! (- 12\|z\|^2+ (G_c(z),p) ) with τ=200τ=200. Comparison We use the pretrained iMF for class-conditional ImageNet generation (Geng et al., 2026b). We compare iMF+SPT+pCN and iMF+SPT+hybrid against their multi-step counterpart SiT+SPT+pCN, where SiT (Ma et al., 2024) is a pretrained class-conditional flow-matching prior integrated using 125 Heun steps. These methods use eight temperature replicas and eight independent ladders. For iMF-based methods, we run 250 adaptation sweeps, 250 burn-in sweeps, and retain 1,000 post-burn-in samples per ladder, yielding 8,000 retained samples per method. SiT+SPT+pCN is extremely slow in this experiment, we evaluate it under a matched per-prompt sampling-time budget equal to the measured hybrid sampling time (more than 1010 hours), with an additional one-hour allowance for initialization, scoring, and image export. We also include best-of-K (Karthik et al., 2023), K=12K=12, as a reward-based selection baseline, retaining the image with the highest CLIP reward among K learned-prior independent samples. We produce 8,0008,000 best-of-K images. Results. Table 2 reports CLIP rewards and sampling costs for two class-prompt pairs, while Figures 2 and 2 show samples (R is the reward) from two ladders per method. JPEG compression is for display only and does not affect reported rewards. The two iMF methods achieve comparable rewards and both outperform SiT+SPT+pCN. Their samples more faithfully depict a running golden retriever amid fallen autumn leaves, whereas SiT+SPT+pCN often misses the autumn setting within its sampling budget of 10-11 hours. For the bald eagle prompt, all three methods depict eagles with outstretched wings, but the iMF methods more consistently include snow-capped mountains whereas SiT+SPT+pCN depicts the eagle against the plain blue sky. In these examples, the hybrid sampler offers no clear improvement in reward or sample quality over iMF+SPT+pCN. Finally, best-of-K achieves the lowest rewards, suggesting guidance is more effective than independent prior draws. The full results are given in Appendix B.3.2. 5 Conclusion We develop few-step iMF source-space inference, linking population training suboptimality and approximation error to posterior stability and end-to-end error bounds under stated assumptions. Hybrid SPT improves efficiency on synthetic benchmarks, while one-step iMF reduces CLIP-guided ImageNet sampling costs relative to multi-step SiT. Quantifying how model-class approximation error ϵapp _app decreases with capacity while controlling the regularity constants remains open. Table 2: CLIP reward and measured sampling cost. Method Mean CLIP reward ↑ Sampling time (h) Mean sweep time (s) Golden retriever: “A photograph of a golden retriever running through fallen autumn leaves.” Best-of-K 0.2864 2.49 – SiT+SPT+pCN 0.3122 10.51 370.89 iMF+SPT+pCN 0.3809 2.40 5.75 iMF+SPT+hybrid 0.3824 10.95 26.28 Bald eagle: “A photograph of a bald eagle soaring above snow-capped mountains, wings spread wide.” Best-of-K 0.3278 2.49 – SiT+SPT+pCN 0.3481 10.70 370.27 iMF+SPT+pCN 0.3936 2.44 5.85 iMF+SPT+hybrid 0.3914 11.15 26.76 (a) SiT + SPT + pCN. (b) iMF + SPT + pCN. (c) iMF + SPT + hybrid. Figure 1: Class: golden retriever, prompt: “A photograph of a golden retriever running through fallen autumn leaves.” (d) SiT + SPT + pCN. (e) iMF + SPT + pCN. (f) iMF + SPT + hybrid. Figure 2: Class: bald eagle, prompt: “A photograph of a bald eagle soaring above snow-capped mountains, wings spread wide.” AI use statement We used generative AI tools (GPT6-Astra, Claude Opus 5, Codex, Copilot Pro+) to assist with language editing; literature search and paper summarization; algorithm implementation, comparative experiments, and code testing; preliminary (unusable) drafts and templates; and developing and checking mathematical arguments; and reviewing the manuscript. We did not use generative AI to generate synthetic datasets. The authors reviewed the AI-assisted material, including text, code, and mathematical arguments. We take full responsibility for the final manuscript, its claims, and all accompanying artifacts, including those produced with generative AI assistance. Acknowledgement This work is supported by NTU-SUG and Singapore Ministry of Education (MOE) AcRF Tier 1 RG17/24 and the Research Council of Finland grant number 369502. We acknowledge computational resources on LUMI supercomputer, owned by the EuroHPC Joint Undertaking and hosted by CSC (Finland) and the LUMI consortium. References Ben-Hamu et al. (2024) H. Ben-Hamu, O. Puny, I. Gat, B. Karrer, U. Singer, and Y. Lipman D-flow: differentiating through flows for controlled generation. arXiv preprint arXiv:2402.14017. Cited by: §1. Beskos et al. (2011) A. Beskos, F. J. Pinski, J. M. Sanz-Serna, and A. M. Stuart Hybrid Monte Carlo on Hilbert spaces. Stochastic Processes and their Applications 121 (10), p. 2201–2230. Cited by: §1, §3.3, §3.3. Beskos et al. (2008) A. Beskos, G. Roberts, A. Stuart, and J. Voss MCMC methods for diffusion bridges. Stochastics and Dynamics 8 (03), p. 319–350. Cited by: §1, §3.3. Betancourt et al. (2014) M. Betancourt, S. Byrne, and M. Girolami Optimizing the integrator step size for Hamiltonian Monte Carlo. arXiv preprint arXiv:1411.6669. Cited by: §B.2. Bissiri et al. (2016) P. G. Bissiri, C. C. Holmes, and S. G. Walker A general framework for updating belief distributions. Journal of the Royal Statistical Society Series B: Statistical Methodology 78 (5), p. 1103–1130. Cited by: §1. Boffi et al. (2024) N. M. Boffi, M. S. Albergo, and E. Vanden-Eijnden Flow map matching with stochastic interpolants: a mathematical framework for consistency models. arXiv preprint arXiv:2406.07507. Cited by: §3.2. Chen et al. (2018) R. T. Chen, Y. Rubanova, J. Bettencourt, and D. K. Duvenaud Neural ordinary differential equations. Advances in neural information processing systems 31. Cited by: §1. Chen et al. (2016) Y. Chen, D. Keyes, K. J. Law, and H. Ltaief Accelerated dimension-independent adaptive metropolis. SIAM Journal on Scientific Computing 38 (5), p. S539–S565. Cited by: §B.2. Cheng et al. (2026) Z. Cheng, L. Wang, and Z. Wang Preconditioned one-step generative modeling for Bayesian inverse problems in function spaces. arXiv preprint arXiv:2603.14798. Cited by: §1. Chung et al. (2023) H. Chung, J. Kim, M. T. Mccann, M. L. Klasky, and J. C. Ye Diffusion posterior sampling for general noisy inverse problems. In The Eleventh International Conference on Learning Representations, Cited by: §B.2, §1, §1, §4.1. Cotter et al. (2013) S. L. Cotter, G. O. Roberts, A. M. Stuart, and D. White MCMC methods for functions: modifying old algorithms to make them faster. Statistical Science, p. 424–446. Cited by: §3.3. Earl and Deem (2005) D. J. Earl and M. W. Deem Parallel tempering: theory, applications, and new perspectives. Physical Chemistry Chemical Physics 7 (23), p. 3910–3916. Cited by: §1. Geng et al. (2026a) Z. Geng, M. Deng, X. Bai, Z. Kolter, and K. He Mean flows for one-step generative modeling. Advances in Neural Information Processing Systems 38, p. 75460–75482. Cited by: §1, §2.2, §3.2. Geng et al. (2026b) Z. Geng, Y. Lu, Z. Wu, E. Shechtman, J. Z. Kolter, and K. He Improved mean flows: on the challenges of fastforward generative models. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, p. 30467–30476. Cited by: §1, §2.2, §3.2, §4.2. Geyer (1991) C. J. Geyer Markov chain Monte Carlo maximum likelihood. Cited by: §3.3. Greenberg et al. (2019) D. Greenberg, M. Nonnenmacher, and J. Macke Automatic posterior transformation for likelihood-free inference. In International conference on machine learning, p. 2404–2414. Cited by: §1. Hairer et al. (1993) E. Hairer, G. Wanner, and S. P. Nørsett Solving ordinary differential equations i: nonstiff problems. Springer. Cited by: §1. Hoffman et al. (2014) M. D. Hoffman A. Gelman et al. The No-U-Turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo. Journal of Machine Learning Research 15 (1), p. 1593–1623. Cited by: §1. Hoffman et al. (2019) M. Hoffman, P. Sountsov, J. V. Dillon, I. Langmore, D. Tran, and S. Vasudevan Neutra-lizing bad geometry in hamiltonian monte carlo using neural transport. arXiv preprint arXiv:1903.03704. Cited by: §1. Karthik et al. (2023) S. Karthik, K. Roth, M. Mancini, and Z. Akata If at first you don’t succeed, try, try again: faithful diffusion-based text-to-image generation by selection. arXiv preprint arXiv:2305.13308. Cited by: §4.2. Korbak et al. (2022) T. Korbak, E. Perez, and C. Buckley RL with KL penalties is better viewed as Bayesian inference. In Findings of the Association for Computational Linguistics: EMNLP 2022, p. 1083–1091. Cited by: §1. Lipman et al. (2023) Y. Lipman, R. T. Q. Chen, H. Ben-Hamu, M. Nickel, and M. Le Flow matching for generative modeling. In The Eleventh International Conference on Learning Representations, Cited by: §1, §2.1, §2.1. Lu and Wang (2026) J. Lu and Z. Wang On variance reduction in learning mean flows. arXiv preprint arXiv:2605.09235. Cited by: footnote 1. Lunz et al. (2018) S. Lunz, O. Öktem, and C. Schönlieb Adversarial regularizers in inverse problems. Advances in neural information processing systems 31. Cited by: §1. Ma et al. (2024) N. Ma, M. Goldstein, M. S. Albergo, N. M. Boffi, E. Vanden-Eijnden, and S. Xie Sit: exploring flow and diffusion-based generative models with scalable interpolant transformers. In European Conference on Computer Vision, p. 23–40. Cited by: §4.2. Otto and Villani (2000) F. Otto and C. Villani Generalization of an inequality by Talagrand and links with the logarithmic Sobolev inequality. Journal of Functional Analysis 173 (2), p. 361–400. Cited by: §A.7. Papamakarios et al. (2021) G. Papamakarios, E. Nalisnick, D. J. Rezende, S. Mohamed, and B. Lakshminarayanan Normalizing flows for probabilistic modeling and inference. Journal of Machine Learning Research 22 (57), p. 1–64. Cited by: §1. Parno and Marzouk (2018) M. D. Parno and Y. M. Marzouk Transport map accelerated markov chain monte carlo. SIAM/ASA Journal on Uncertainty Quantification 6 (2), p. 645–682. Cited by: §1. Patel and Oberai (2019) D. Patel and A. A. Oberai Bayesian inference with generative adversarial network priors. arXiv preprint arXiv:1907.09987. Cited by: §1. Patel et al. (2022) D. V. Patel, D. Ray, and A. A. Oberai Solution of physics-based bayesian inverse problems with deep generative priors. Computer Methods in Applied Mechanics and Engineering 400, p. 115428. Cited by: §1. Purohit et al. (2025) V. Purohit, M. Repasky, J. Lu, Q. Qiu, Y. Xie, and X. Cheng Consistency posterior sampling for diverse image synthesis. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), p. 28327–28336. Cited by: §1, §3.1. Radev et al. (2023) S. T. Radev, M. Schmitt, L. Schumacher, L. Elsemüller, V. Pratz, Y. Schälte, U. Köthe, and P. Bürkner BayesFlow: amortized bayesian workflows with neural networks. arXiv preprint arXiv:2306.16015. Cited by: §1. Radford et al. (2021) A. Radford, J. W. Kim, C. Hallacy, A. Ramesh, G. Goh, S. Agarwal, G. Sastry, A. Askell, P. Mishkin, J. Clark, et al. Learning transferable visual models from natural language supervision. In International conference on machine learning, p. 8748–8763. Cited by: §1, §4.2. Ramzi et al. (2020) Z. Ramzi, B. Remy, F. Lanusse, J. Starck, and P. Ciuciu Denoising score-matching for uncertainty quantification in inverse problems. arXiv preprint arXiv:2011.08698. Cited by: §1. Stuart (2010) A. M. Stuart Inverse problems: a Bayesian perspective. Acta numerica 19, p. 451–559. Cited by: §1. Tong et al. (2024) A. Tong, K. Fatras, N. Malkin, G. Huguet, Y. Zhang, J. Rector-Brooks, G. Wolf, and Y. Bengio Improving and generalizing flow-based generative models with minibatch optimal transport. Transactions on Machine Learning Research. Cited by: §1, §1. Villani et al. (2009) C. Villani et al. Optimal transport: old and new. Vol. 338, Springer. Cited by: §3.1. Vincent (2011) P. Vincent A connection between score matching and denoising autoencoders. Neural computation 23 (7), p. 1661–1674. Cited by: §1, §1. Wang et al. (2026a) S. Wang, J. Keller, T. Transue, D. B. Brown, T. Strohmer, and B. Wang Test-time guidance for flow-based generative models via parallel tempering on source distributions. In International Conference on Machine Learning, Cited by: §A.7, §A.7, §1, §1, §1, §3.3, §3.3, §3.3. Wang et al. (2026b) Z. Wang, A. Harting, M. Barreau, M. Zavlanos, and K. H. Johansson Source-guided flow matching. In International Conference on Learning Representations, Vol. 2026, p. 106484–106514. Cited by: §1, §1, §3.1. Zhang et al. (2025) B. J. Zhang, Y. M. Marzouk, and K. Spiliopoulos Transport map unadjusted Langevin algorithms: learning and discretizing perturbed samplers. Foundations of Data Science 7 (3), p. 705–736. Cited by: §1. Appendix A Proofs and clarifications A.1 Equivalence between ℒiMFL_iMF and ℛiMFR_iMF Write ℒiMF(θ)=‖Aθ−V‖22L_iMF(θ)=E\|A_θ-V\|_2^2, where the predictor Aθ:=uθ(WT,R,T)+(T−R)[∂tuθ(WT,R,T)+(∂wuθ)(WT,R,T)vθ(WT,T)]A_θ:=u_θ(W_T,R,T)+(T-R) [ _tu_θ(W_T,R,T)+( _wu_θ)(W_T,R,T)\,v_θ(W_T,T) ] is measurable with respect to (WT,R,T)(W_T,R,T). Since (R,T)(R,T) is independent of (X,Z)(X,Z), [V∣WT,R,T]=v(WT,T).E[V W_T,R,T]=v(W_T,T). Hence ⟨Aθ−v(WT,T),V−v(WT,T)⟩=0.E A_θ-v(W_T,T),\,V-v(W_T,T) =0. Expanding the squared norm therefore gives ℒiMF(θ)=ℛiMF(θ)+‖V−v(WT,T)‖22⏟C,L_iMF(θ)=R_iMF(θ)+ E\|V-v(W_T,T)\|_2^2_C, where C is independent of θ. A.2 Proof of relation (14) Since ℛ(θ)=ℒ(θ)+CR(θ)=L(θ)+C, it holds infℛ=infℒ+C = +C. Therefore, ℛ(θ)=ℒ(θ)−infℒ+infℛR(θ)=L(θ)- + . A.3 Details on Remark 3.4 Set A:=sup0≤r≤t≤1‖q¯θ(⋅,r,t)‖L2(pt),B:=sup0≤r≤t≤1‖∂rq¯θ(⋅,r,t)‖L2(pt).A:= _0≤ r≤ t≤ 1\| q_θ(·,r,t)\|_L^2(p_t), B:= _0≤ r≤ t≤ 1\| _r q_θ(·,r,t)\|_L^2(p_t). By hypothesis, A,B<∞A,B<∞. For 0≤r<s≤t≤10≤ r<s≤ t≤ 1, the fundamental theorem of calculus and Minkowski’s inequality give ‖q¯θ(⋅,r,t)−q¯θ(⋅,s,t)‖L2(pt)≤∫rs‖∂rq¯θ(⋅,a,t)‖L2(pt)a≤B(s−r).\| q_θ(·,r,t)- q_θ(·,s,t)\|_L^2(p_t)≤ _r^s\| _r q_θ(·,a,t)\|_L^2(p_t)\,da≤ B(s-r). Consequently, |‖q¯θ(⋅,r,t)‖L2(pt)2−‖q¯θ(⋅,s,t)‖L2(pt)2|≤2AB(s−r). |\| q_θ(·,r,t)\|_L^2(p_t)^2-\| q_θ(·,s,t)\|_L^2(p_t)^2 |≤ 2AB(s-r). Splitting the integrals defining Gθ(r)G_θ(r) and Gθ(s)G_θ(s) therefore yields |Gθ(r)−Gθ(s)| |G_θ(r)-G_θ(s)| ≤∫rs‖q¯θ(⋅,r,t)‖L2(pt)2t ≤ _r^s\| q_θ(·,r,t)\|_L^2(p_t)^2\,dt +∫s1|∥q¯θ(⋅,r,t)∥L2(pt)2−∥q¯θ(⋅,s,t)∥L2(pt)2|dt + _s^1 |\| q_θ(·,r,t)\|_L^2(p_t)^2-\| q_θ(·,s,t)\|_L^2(p_t)^2 |\,dt ≤A2(s−r)+2AB(1−s)(s−r) ≤ A^2(s-r)+2AB(1-s)(s-r) ≤(A2+2AB)|s−r|. ≤(A^2+2AB)|s-r|. Now suppose that both uθu_θ and vθv_θ are finite compositions of affine maps and coordinatewise SiLU activations, with direct spatial and time inputs. The SiLU activation ρ(a)=a/(1+e−a)ρ(a)=a/(1+e^-a) is smooth and has bounded first and second derivatives. By the chain rule, these networks have bounded first and second input derivatives and at most linear growth. Since the time inputs lie in [0,1][0,1], there is a finite constant CθC_θ such that ‖uθ(w,r,t)‖+‖vθ(w,t)‖≤Cθ(1+‖w‖),\|u_θ(w,r,t)\|+\|v_θ(w,t)\|≤ C_θ(1+\|w\|), and the first and second input derivatives of uθu_θ are uniformly bounded. Because v and vθv_θ do not depend on r, ∂rq¯θ=∂ruθ−∂tuθ−(∂wuθ)vθ+(t−r)[∂rtuθ+(∂rwuθ)vθ]. _r q_θ= _ru_θ- _tu_θ-( _wu_θ)v_θ+(t-r) [ _rtu_θ+( _rwu_θ)v_θ ]. The preceding network bounds imply, after enlarging CθC_θ if necessary, ‖q¯θ(w,r,t)‖ \| q_θ(w,r,t)\| ≤Cθ(1+‖w‖)+‖v(w,t)‖, ≤ C_θ(1+\|w\|)+\|v(w,t)\|, ‖∂rq¯θ(w,r,t)‖ \| _r q_θ(w,r,t)\| ≤Cθ(1+‖w‖), ≤ C_θ(1+\|w\|), uniformly over 0≤r≤t≤10≤ r≤ t≤ 1. Let M2:=‖X‖2<∞M_2:=E\|X\|^2<∞ where X∼π0X _0 (Assumption 3.1). Independence of X and Z∼(0,Id)Z (0,I_d) gives ‖Wt‖2=(1−t)2M2+t2d≤M2+d.E\|W_t\|^2=(1-t)^2M_2+t^2d≤ M_2+d. Moreover, since v(Wt,t)=[Z−X∣Wt]v(W_t,t)=E[Z-X W_t], conditional Jensen’s inequality gives ‖v(⋅,t)‖L2(pt)2≤‖Z−X‖2=M2+d.\|v(·,t)\|_L^2(p_t)^2 \|Z-X\|^2=M_2+d. These estimates establish A,B<∞A,B<∞, proving the network claim. The resulting constants may depend on the network parameters. A.4 Details on Remark 3.6 Write eθ(w,t):=v(w,t)−vθ(w,t)e_θ(w,t):=v(w,t)-v_θ(w,t). The operator-norm bound implies, for r<tr<t, (t−r)‖∂wuθ(w,r,t)eθ(w,t)‖≤Bu‖eθ(w,t)‖.(t-r)\| _wu_θ(w,r,t)e_θ(w,t)\|≤ B_u\|e_θ(w,t)\|. Squaring, integrating against ptp_t, and summing over the schedule intervals gives ∑k=1K∫tk−1tk(t−tk−1)2‖∂wuθ(⋅,tk−1,t)eθ(⋅,t)‖L2(pt)2t _k=1^K _t_k-1^t_k(t-t_k-1)^2\| _wu_θ(·,t_k-1,t)e_θ(·,t)\|_L^2(p_t)^2\,dt ≤Bu2∑k=1K∫tk−1tk‖eθ(⋅,t)‖L2(pt)2t=Bu2ℛFM(θ), ≤ B_u^2 _k=1^K _t_k-1^t_k\|e_θ(·,t)\|_L^2(p_t)^2\,dt=B_u^2R_FM(θ), where the last equality uses T∼Unif[0,1]T [0,1] and the fact that the schedule intervals partition [0,1][0,1]. If uθu_θ is uniformly LuL_u-Lipschitz in w, then ‖∂wuθ(w,r,t)‖op≤Lu\| _wu_θ(w,r,t)\|_op≤ L_u. Since t−r≤1t-r≤ 1, the sufficient condition follows with Bu=LuB_u=L_u. Finally, (15) controls only the Jacobian acting on eθe_θ, averaged over space and time; it does not require a uniform operator-norm bound. For example, when vθ=v_θ=v, its left-hand side vanishes. For finite affine-SiLU networks with direct spatial and time inputs, the spatial Jacobian is uniformly bounded (similar to the argument in Appendix A.3): Lu:=supw∈ℝd,(r,t)∈[0,1]2‖∂wuθ(w,r,t)‖op<∞.L_u:= _w ^d,\,(r,t)∈[0,1]^2\| _wu_θ(w,r,t)\|_op<∞. A.5 Proof of Theorem 3.7 Population risk control: By (14) and the definition of the joint marginal risk, ℛiMF(θ)+λℛFM(θ)=ℛ(θ)=δ.R_iMF(θ)+ _FM(θ)=R(θ)=δ. Therefore, ℛiMF(θ)≤δR_iMF(θ)≤δ and ℛFM(θ)≤δ/λR_FM(θ)≤δ/λ. Recovering fixed-time slices: Independent uniform time sampling on [0,1]2[0,1]^2 implies ∫01Gθ(r)r=∫01∫r1‖q¯θ(⋅,r,t)‖L2(pt)2tr≤∫01∫01‖q¯θ(⋅,r,t)‖L2(pt)2tr=ℛiMF(θ). _0^1G_θ(r)\,dr= _0^1 _r^1\| q_θ(·,r,t)\|_L^2(p_t)^2\,dt\,dr≤ _0^1 _0^1\| q_θ(·,r,t)\|_L^2(p_t)^2\,dt\,dr=R_iMF(θ). Fix r∈[0,1]r∈[0,1] and h∈(0,1]h∈(0,1]. Choose an interval I=[a,b]⊂[0,1]I=[a,b]⊂[0,1] such that a≤r≤ba≤ r≤ b and b−a=hb-a=h. By Assumption 3.3, Gθ(r)≤Gθ(s)+LG|s−r|,s∈I.G_θ(r)≤ G_θ(s)+L_G|s-r|, s∈ I. Averaging this inequality over I gives Gθ(r)≤1h∫IGθ(s)s+LGh∫I|s−r|s.G_θ(r)≤ 1h _IG_θ(s)\,ds+ L_Gh _I|s-r|\,ds. The first term is at most δ/hδ/h, since Gθ≥0G_θ≥ 0 and ∫01Gθ(s)s≤ℛiMF(θ)≤δ _0^1G_θ(s)\,ds _iMF(θ)≤δ. For the second term, note that ∫ab|s−r|s=(r−a)2+(b−r)22≤((r−a)+(b−r))22=h22. _a^b|s-r|\,ds= (r-a)^2+(b-r)^22≤ ((r-a)+(b-r))^22= h^22. Consequently, Gθ(r)≤δh+LGh2. G_θ(r)≤ δh+ L_Gh2. (22) If 0<δ≤LG/20<δ≤ L_G/2, choosing h=2δ/LG∈(0,1]h= 2δ/L_G∈(0,1] to minimize the bound in (22), we get Gθ(r)≤(2LGδ)1/2,∀r∈[0,1].G_θ(r)≤(2L_Gδ)^1/2, ∀ r∈[0,1]. (23) If δ=0δ=0, the same bound follows by letting h↓0h 0 in (22). Local transport error: Define qθ(w,r,t) q_θ(w,r,t) :=uθ(w,r,t)−v(w,t)+(t−r)[∂tuθ(w,r,t)+∂wuθ(w,r,t)v(w,t)] :=u_θ(w,r,t)-v(w,t)+(t-r) [ _tu_θ(w,r,t)+ _wu_θ(w,r,t)v(w,t) ] =q¯θ(w,r,t)+(t−r)∂wuθ(w,r,t)(v(w,t)−vθ(w,t)). = q_θ(w,r,t)+(t-r) _wu_θ(w,r,t) (v(w,t)-v_θ(w,t) ). Fix 0≤r<t≤10≤ r<t≤ 1 and, for ptp_t-almost every w, set (see Assumption 3.1 for the flow definition) ws:=φs,t(w),Ds:=(s−r)uθ(ws,r,s)−ws,s∈[r,t].w_s:= _s,t(w), D_s:=(s-r)u_θ(w_s,r,s)-w_s, s∈[r,t]. Assumptions 3.1 and 3.2 ensure that DsD_s is absolutely continuous. The chain rule gives dsDs=qθ(ws,r,s) ddsD_s=q_θ(w_s,r,s) for almost every s. Since Dr=−φr,t(w)D_r=- _r,t(w) and Dt=(t−r)uθ(w,r,t)−wD_t=(t-r)u_θ(w,r,t)-w, integration gives (t−r)uθ(w,r,t)−w+φr,t(w)=∫rtqθ(φs,t(w),r,s)s.(t-r)u_θ(w,r,t)-w+ _r,t(w)= _r^tq_θ( _s,t(w),r,s)\,ds. (24) For the fixed generation schedule, let Fk:=φtk−1,tkF_k:= _t_k-1,t_k be the exact step map, and define Qk Q_k :=∫tk−1tk‖q¯θ(⋅,tk−1,t)‖L2(pt)2t, := _t_k-1^t_k\| q_θ(·,t_k-1,t)\|_L^2(p_t)^2\,dt, Jv,k J_v,k :=∫tk−1tk(t−tk−1)2‖∂wuθ(⋅,tk−1,t)(v−vθ)(⋅,t)‖L2(pt)2t. := _t_k-1^t_k(t-t_k-1)^2\| _wu_θ(·,t_k-1,t)(v-v_θ)(·,t)\|_L^2(p_t)^2\,dt. Equation (24) gives Fk,θ(w)−Fk(w)=−∫tk−1tkqθ(φs,tk(w),tk−1,s)ds.F_k,θ(w)-F_k(w)=- _t_k-1^t_kq_θ( _s,t_k(w),t_k-1,s)\,ds. Using (φs,tk)#ptk=ps( _s,t_k)_\#p_t_k=p_s, Minkowski’s inequality, and Cauchy-Schwarz in time, we obtain ‖Fk,θ−Fk‖L2(ptk) \|F_k,θ-F_k\|_L^2(p_t_k) ≤∫tk−1tk‖qθ(φs,tk(⋅),tk−1,s)‖L2(ptk)s ≤ _t_k-1^t_k\|q_θ( _s,t_k(·),t_k-1,s)\|_L^2(p_t_k)\,ds =∫tk−1tk‖qθ(⋅,tk−1,s)‖L2(ps)s = _t_k-1^t_k\|q_θ(·,t_k-1,s)\|_L^2(p_s)\,ds ≤∫tk−1tk‖q¯θ(⋅,tk−1,s)‖L2(ps)s ≤ _t_k-1^t_k\| q_θ(·,t_k-1,s)\|_L^2(p_s)\,ds +∫tk−1tk(s−tk−1)∥∂wuθ(⋅,tk−1,s)(v−vθ)(⋅,s)∥L2(ps)ds + _t_k-1^t_k(s-t_k-1)\| _wu_θ(·,t_k-1,s)(v-v_θ)(·,s)\|_L^2(p_s)\,ds ≤hk(Qk+Jv,k). ≤ h_k ( Q_k+ J_v,k ). (25) Composing the step maps: The flow property gives S=F1∘⋯∘FK=φ0,1S=F_1 ·s F_K= _0,1. Define the intermediate compositions Hk:=F1,θ∘⋯∘Fk,θ∘Fk+1∘⋯∘FK,k=0,…,K,H_k:=F_1,θ ·s F_k,θ F_k+1 ·s F_K, k=0,…,K, where empty compositions are the identity. Thus H0=SH_0=S, HK=SθH_K=S_θ, and Sθ(z)−S(z)=∑k=1K(Hk(z)−Hk−1(z)).S_θ(z)-S(z)= _k=1^K (H_k(z)-H_k-1(z) ). The maps HkH_k and Hk−1H_k-1 differ only at the kkth step. Their common inner composition is the exact flow Fk+1∘⋯∘FK=φtk,1F_k+1 ·s F_K= _t_k,1, while their common outer composition F1,θ∘⋯∘Fk−1,θF_1,θ ·s F_k-1,θ has Lipschitz constant at most Λk=∏j<kLj _k= _j<kL_j. Hence ‖Hk(z)−Hk−1(z)‖≤Λk‖Fk,θ(φtk,1(z))−Fk(φtk,1(z))‖.\|H_k(z)-H_k-1(z)\|≤ _k\|F_k,θ( _t_k,1(z))-F_k( _t_k,1(z))\|. Since (φtk,1)#γ=ptk( _t_k,1)_\#γ=p_t_k, taking L2(γ)L^2(γ) norms and applying the triangle inequality yields ‖Sθ−S‖L2(γ)≤∑k=1KΛk‖Fk,θ−Fk‖L2(ptk).\|S_θ-S\|_L^2(γ)≤ _k=1^K _k\|F_k,θ-F_k\|_L^2(p_t_k). Consequently, ‖Sθ−S‖L2(γ) \|S_θ-S\|_L^2(γ) ≤∑k=1KΛk‖Fk,θ−Fk‖L2(ptk) ≤ _k=1^K _k\|F_k,θ-F_k\|_L^2(p_t_k) ≤∑k=1KΛkhk(Qk+Jv,k) ≤ _k=1^K _k h_k ( Q_k+ J_v,k ) ≤AK(2LGδ)1/4+BK(∑k=1KJv,k)1/2 ≤ A_K(2L_Gδ)^1/4+B_K ( _k=1^KJ_v,k )^1/2 ≤AK(2LGδ)1/4+BuBKλδ1/2. ≤ A_K(2L_Gδ)^1/4+ B_uB_K λδ^1/2. The third line uses Cauchy–Schwarz over k for the Jv,kJ_v,k terms. Squaring proves (16). For K=1K=1, we have h1=Λ1=A1=B1=1h_1= _1=A_1=B_1=1, so the same argument gives the one-step transport bound. In this case there is no outer learned map in the composition estimate, and no Lipschitz condition on F1,θF_1,θ is needed. Finally, if uθu_θ is uniformly LuL_u-Lipschitz in space, then Lk≤1+hkLu,Λk≤exp(Lu∑j<khj)≤eLu.L_k≤ 1+h_kL_u, _k≤ (L_u _j<kh_j )≤ e^L_u. Since ∑khk=1 _kh_k=1, Cauchy–Schwarz yields AK≤eLuKA_K≤ e^L_u K and BK≤eLuB_K≤ e^L_u. Also LSθ≤∏kLk≤eLuL_S_θ≤ _kL_k≤ e^L_u, as used in the data-space posterior bounds. All estimates concern the fixed generation schedule. A.6 Proof of Theorem 3.9 Write e:=ℰ1/2e:=E^1/2, a(z):=ℓy(Sθ(z))a(z):= _y(S_θ(z)) and b(z):=ℓy(S(z))b(z):= _y(S(z)). The weight assumptions give |Zθ,y−Zy|≤∫|a−b|γ≤Lye.|Z_θ,y-Z_y|≤ |a-b|\,dγ≤ L_ye. Hence Lye≤Zy/2L_ye≤ Z_y/2 ensures Zθ,y≥Zy/2>0Z_θ,y≥ Z_y/2>0. Let U2(z)=1+‖z‖2U_2(z)=1+\|z\|^2 and use the weighted-variation convention ‖μ−ν‖U2:=sup|f|≤U2|∫fμ−∫fν|.\|μ-ν\|_U_2:= _|f|≤ U_2 | f\,dμ- f\,dν |. For each such f, decompose the normalized densities as ∫fd(νθy−νy)=1Zθ,y∫f(a−b)γ+Zy−Zθ,yZyZθ,y∫fbγ. f\,d( _θ^y-ν^y)= 1Z_θ,y f(a-b)\,dγ+ Z_y-Z_θ,yZ_yZ_θ,y fb\,dγ. Since γU2=1+dE_γU_2=1+d and γU22=d2+4d+1E_γU_2^2=d^2+4d+1, Cauchy-Schwarz gives ‖νθy−νy‖U2 \| _θ^y-ν^y\|_U_2 ≤LyeZθ,y[d2+4d+1+My(1+d)Zy] ≤ L_yeZ_θ,y [ d^2+4d+1+ M_y(1+d)Z_y ] ≤Cy(d)2e. ≤ C_y(d)2e. For completeness, couple the common part of two measures identically and couple their remaining positive parts. The inequality ‖z−z′‖2≤2‖z‖2+2‖z′‖2\|z-z \|^2≤ 2\|z\|^2+2\|z \|^2 then gives W22(μ,ν)≤2∫‖z‖2d|μ−ν|(z)≤2‖μ−ν‖U2.W_2^2(μ,ν)≤ 2 \|z\|^2\,d|μ-ν|(z)≤ 2\|μ-ν\|_U_2. This proves (18). All source second moments are finite because the normalized weights are bounded relative to γ. For the data-space statement, first note that S#γ=π0S_\#γ= _0 implies S#νy=πyS_\#ν^y=π^y by the definition of a pushforward. Insert the intermediate measure (Sθ)#νy(S_θ)_\#ν^y. Lipschitz propagation and the triangle inequality give W2((Sθ)#νθy,πy)≤LSθW2(νθy,νy)+W2((Sθ)#νy,S#νy).W_2 ((S_θ)_\# _θ^y,π^y )≤ L_S_θW_2( _θ^y,ν^y)+W_2 ((S_θ)_\#ν^y,S_\#ν^y ). The common-input coupling bounds the second term by W22((Sθ)#νy,S#νy)≤νy‖Sθ−S‖2≤MyZyℰ,W_2^2 ((S_θ)_\#ν^y,S_\#ν^y ) _ν^y\|S_θ-S\|^2≤ M_yZ_yE, which proves the data-space stability. Substituting (17) yields (19)–(20). A.7 Geometric convergence of the source samplers Fix θ and y, and write ω(z):=ℓy(Sθ(z)),μj(z):=Zj−1ω(z)βjγ(z),Zj:=∫ωβjγ,Π:=⨂j=0Jμj.ω(z):= _y(S_θ(z)), _j(dz):=Z_j^-1ω(z) _jγ(dz), Z_j:= ω _j\,dγ, := _j=0^J _j. Thus μ0=γ _0=γ and μJ=νθy _J= _θ^y. Let ηn _n be the joint law after n production sweeps and νn _n its cold marginal. We use χ2(η∥Π):=∫(dη/dΠ−1)2dΠχ^2(η\| ):= (dη/d -1)^2\,d . Assumption A.1 (Source weight and production sampler). The source weight is measurable and satisfies 0<m≤ω(z)≤M<∞0<m≤ω(z)≤ M<∞ for every z∈ℝdz ^d, where m<Mm<M. The ladder 0=β0<⋯<βJ=10= _0<·s< _J=1, J≥1J≥ 1, and all local-kernel parameters are fixed during production. Each sweep independently applies the local kernels to all replicas, then attempts each adjacent Metropolis swap once in a prescribed order; this order may vary with the sweep index. Every pCN angle satisfies hβj∈(0,π/2]h_ _j∈(0,π/2]. For the hybrid sampler, Ψθ=−logω _θ=- ω is continuously differentiable, and the cold update is the Metropolis-corrected split HMC kernel of Section 3.3, with fresh independent Gaussian momentum at each update. The initial joint law satisfies D0:=χ2(η0∥Π)<∞D_0:=χ^2( _0\| )<∞. The positive lower bound is an additional condition for the convergence result. It holds, for example, for exponentiated bounded rewards. Geometric convergence of SPT+pCN is established by Wang et al. (2026a). The following statement gives explicit bounds for both local-kernel choices and extends the guarantee to the hybrid sampler. Theorem A.2 (Geometric convergence of SPT+pCN and SPT+hybrid). Under Assumption A.1, set κ:=M/mκ:=M/m. For iMF+SPT+pCN, define gpCN:=min0≤j≤Jκ−βj(1−coshβj),qpCN:=1−gpCN∈(0,1).g_pCN:= _0≤ j≤ Jκ^- _j(1- h_ _j), q_pCN:=1-g_pCN∈(0,1). Then, for every integer n≥0n≥ 0, χ2(ηn∥Π)≤D0qpCN2n,W2(νn,νθy)≤2κD0qpCNn.χ^2( _n\| )≤ D_0q_pCN^2n, W_2( _n, _θ^y)≤ 2κ D_0\,q_pCN^n. (26) For iMF+SPT+hybrid, define ghot:=min0≤j<Jκ−βj(1−coshβj),a:=κ−(1−βJ−1),qhyb:=1−ghot2a∈(0,1).g_hot:= _0≤ j<Jκ^- _j(1- h_ _j), a:=κ^-(1- _J-1), q_hyb:=1-g_hot^2a∈(0,1). Then χ2(ηn∥Π) χ^2( _n\| ) ≤D0qhyb2⌊n/2⌋, ≤ D_0q_hyb^2 n/2 , (27) W2(νn,νθy) W_2( _n, _θ^y) ≤2κD0qhyb⌊n/2⌋. ≤ 2κ D_0\,q_hyb n/2 . The hybrid bounds remain valid if the cold update is replaced by any fixed μJ _J-invariant Markov kernel. In particular, no separate spectral-gap assumption on split HMC is required. If η0=γ⊗(J+1) _0=γ (J+1), then D0≤κ∑j=0Jβj−1D_0≤κ _j=0^J _j-1. Proof. All operator norms below are on the centered space L02(Π):=f∈L2(Π):∫fΠ=0L_0^2( ):=\f∈ L^2( ): f\,d =0\, unless another measure is specified. Each local kernel and each Metropolis swap preserves Π , so each is an L2(Π)L^2( ) contraction. Their composition therefore preserves the joint target, whose cold marginal is νθy _θ^y. Contraction of a pCN update. Let j Q_j be the Gaussian pCN proposal kernel and j A_j its Metropolis correction for μj _j. The Gaussian kernel is self-adjoint and positive on L2(γ)L^2(γ), with centered norm coshβj h_ _j. Writing aj:=dμj/dγa_j:=d _j/dγ, its accepted transition measure is minaj(z),aj(z′)γ(dz)j(z,dz′). \a_j(z),a_j(z )\\,γ(dz) Q_j(z,dz ). Hence comparison of Dirichlet forms and variances gives ⟨f,(I−j)f⟩μj f,(I- A_j)f _ _j ≥(infaj)⟨f,(I−j)f⟩γ ≥( a_j) f,(I- Q_j)f _γ ≥infajsupaj(1−coshβj)Varμj(f) ≥ a_j a_j(1- h_ _j)Var_ _j(f) ≥κ−βj(1−coshβj)Varμj(f). ≥κ^- _j(1- h_ _j)Var_ _j(f). The corrected kernel is also positive. Indeed, the identity minaj(z),aj(z′)=∫0∞aj(z)≥saj(z′)≥sds \a_j(z),a_j(z )\= _0^∞1_\a_j(z)≥ s\1_\a_j(z )≥ s\\,ds expresses the accepted part of its quadratic form as an integral of nonnegative quadratic forms of j Q_j; the rejection part is nonnegative as well. Reversibility and positivity therefore imply ‖j‖L02(μj)→L02(μj)≤1−κ−βj(1−coshβj).\| A_j\|_L_0^2( _j)→ L_0^2( _j)≤ 1-κ^- _j(1- h_ _j). (28) All-pCN sweeps. When every replica uses pCN, tensorization of (28) bounds the centered norm of the joint local update by qpCNq_pCN. The subsequent swap sweep is a contraction, so the centered norm after n sweeps is at most qpCNnq_pCN^n. Applying the adjoint operator to dη0/dΠ−1d _0/d -1 proves the first bound in (26). Hybrid sweeps. Let H be the cold invariant kernel and write =(⨂j<Jj)⊗ L=( _j<J A_j) H for the local update. Let E be the orthogonal projection onto functions of the cold coordinate, obtained by integrating all hot coordinates against ⨂j<Jμj _j<J _j. The operators L and E commute, and (28) gives ‖≤1,‖(I−)‖≤r,r:=1−ghot∈[0,1).\| L E\|≤ 1, \| L(I- E)\|≤ r, r:=1-g_hot∈[0,1). For a complete swap sweep S, every swap involving only hot replicas satisfies j=j= E S_j= S_j E= E. Since the edge (J−1,J)(J-1,J) is attempted once, =J−1 E S E= E S_J-1 E. On functions of the cold coordinate this is an independence Metropolis kernel with target μJ _J and proposal μJ−1 _J-1. Moreover, dμJdμJ−1(z)=ZJ−1ZJω(z)1−βJ−1≤κ1−βJ−1=a−1. d _Jd _J-1(z)= Z_J-1Z_Jω(z)^1- _J-1≤κ^1- _J-1=a^-1. Its transition kernel thus minorizes aμJa _J, yielding ‖≤1−a\| E S E\|≤ 1-a on L02(Π)L_0^2( ). Set :=+r(I−) D:= E+r(I- E). For r>0r>0, the operator :=+r−1(I−) T:= L E+r^-1 L(I- E) is a contraction and satisfies == L= T D= D T. For r=0r=0 the same factorization holds with = T= L. Consequently, ‖ \| L S L\| ≤‖ ≤\| D S D\| ≤r2+2r(1−r)+(1−r)2(1−a) ≤ r^2+2r(1-r)+(1-r)^2(1-a) =1−ghot2a=qhyb. =1-g_hot^2a=q_hyb. Here we expanded =rI+(1−r) D=rI+(1-r) E and used that S, E S, and S E are contractions. Each production sweep acts on functions as n L S_n, so every two consecutive sweeps contract the centered norm by at most qhybq_hyb, regardless of their prescribed swap orders. Applying the adjoint product to the initial density difference proves the χ2χ^2 bound in (27). This argument only uses invariance of H, which proves the assertion for general cold kernels. Cold-marginal Wasserstein bounds. The density of μJ _J relative to γ has upper-to-lower ratio at most κ. The Gaussian logarithmic Sobolev inequality, bounded perturbation, and the implication from logarithmic Sobolev to transport inequalities give (Otto and Villani, 2000) W22(ν,μJ)≤2κKL(ν∥μJ)≤2κχ2(ν∥μJ);W_2^2(ν, _J)≤ 2 (ν\| _J)≤ 2κχ^2(ν\| _J); see also Wang et al. (2026a, Lemma A.10). Conditional Jensen’s inequality gives χ2(νn∥μJ)≤χ2(ηn∥Π)χ^2( _n\| _J)≤χ^2( _n\| ). Combining these inequalities with the joint bounds proves (26) and (27). In particular, the latter has the form CρnCρ^n with ρ=qhybρ= q_hyb and C=2κD0/qhybC= 2κ D_0/q_hyb. Finally, for independent Gaussian initialization, D0+1=∏j=0J[(∫ωβjγ)(∫ω−βjγ)]≤κ∑j=0Jβj.D_0+1= _j=0^J [ ( ω _j\,dγ ) ( ω^- _j\,dγ ) ]≤κ _j=0^J _j. ∎ Corollary A.3 (End-to-end convergence). Suppose the assumptions of Theorem 3.7, Proposition 3.9, and Theorem A.2 hold. Assume also that SθS_θ is LSθL_S_θ-Lipschitz, 0≤δ≤LG/20≤δ≤ L_G/2, and LybK(δ)≤Zy/2L_yb_K(δ)≤ Z_y/2. Define sn:=2κD0qpCNn,for iMF+SPT+pCN,qhyb⌊n/2⌋,for iMF+SPT+hybrid.s_n:= 2κ D_0 casesq_pCN^n,&for iMF+SPT+pCN,\\ q_hyb n/2 ,&for iMF+SPT+hybrid. cases Then the data-space output law satisfies W2((Sθ)#νn,πy)≤ W_2 ((S_θ)_\# _n,π^y )≤ LSθsn+LSθCy(d)1/2bK(δ)1/2 L_S_θs_n+L_S_θC_y(d)^1/2b_K(δ)^1/2 (29) +MyZybK(δ). + M_yZ_yb_K(δ). Appendix B Additional experiments B.1 Empirical loss-to-map-error relationship We examine the relationship between the marginal joint risk and transport-map error in Theorem 3.7 using the two-dimensional banana prior X=(U, 0.62(U2−1)+0.30V),U,V∼iid(0,1).X= (U,\;0.62(U^2-1)+0.30V ), U,V iid N(0,1). For independent X∼π0X _0 and Z∼(0,I2)Z (0,I_2), we use the interpolation Wt=(1−t)X+tZW_t=(1-t)X+tZ and marginal velocity v(w,t)=[Z−X∣Wt=w]v(w,t)=E[Z-X W_t=w]. The reference map S=φ0,1S= _0,1 is the transport induced by this marginal velocity, rather than the explicit triangular map used to generate banana samples. Numerical reference. Conditioning on the first coordinate of WtW_t gives a Gaussian conditional distribution for U. Integrating out V analytically reduces evaluation of v(w,t)v(w,t) to one-dimensional Gaussian quadrature, with weights determined by the second coordinate of WtW_t. We use 128-point Gauss–Hermite quadrature and integrate the marginal ODE from t=1t=1 to t=0t=0 using DOP853 in double precision, with relative and absolute tolerances 10−810^-8 and 10−1010^-10. Doubling the quadrature order changes the reference velocities by 6.9×10−116.9× 10^-11 RMS over the risk-evaluation inputs. On 128 source samples, doubling the quadrature order and tightening both ODE tolerances tenfold change the reference maps by 8.3×10−158.3× 10^-15 and 8.7×10−98.7× 10^-9 RMS, respectively. These sensitivity checks indicate that numerical reference error is negligible at the observed map-error scale. Training and evaluation. We train an iMF model with three hidden SiLU layers of width 128 and separate average-velocity and auxiliary-velocity heads on 20,000 fixed prior samples. Training uses Adam with learning rate 10−310^-3, weight decay 10−510^-5, batch size 512, and 16,000 updates. We use independent data–noise coupling, independent unsorted uniform times (r,t)∈[0,1]2(r,t)∈[0,1]^2, and auxiliary-loss weight λ=1λ=1. The figure reports a single training run with seed 0. At initialization and every 500 updates, we estimate δ=RiMF+λRFMδ=R_iMF+λ R_FM directly from the marginal residuals using 8,192 held-out interpolation samples and the numerical reference velocity. Thus, δ δ estimates the marginal risk appearing in Theorem 3.7, rather than the raw training loss containing irreducible regression noise. We also estimate the squared map error E^K=12048∑i=12048‖Sθ,K(zi)−S(zi)‖22,zi∼(0,I2), E_K= 12048 _i=1^2048 \|S_θ,K(z_i)-S(z_i) \|_2^2, z_i (0,I_2), for K∈1,6K∈\1,6\ uniformly spaced generation steps. The evaluation draws are independent of the training dataset and fixed across checkpoints. Results. Figure 3 shows that decreasing marginal risk is generally accompanied by decreasing map error for both generation schedules. Over training, δ δ decreases from 2.7432.743 to 0.06530.0653, while the squared map error decreases from 0.5540.554 to 0.01980.0198 for one-step generation and from 0.5540.554 to 0.007320.00732 for six-step generation. The relationship is not strictly monotonic across checkpoints. These results provide empirical support for the loss-to-map-error connection in Theorem 3.7, but do not establish its asymptotic rate or uniform control of the regularity constants, which may change during training. Figure 3: Marginal joint risk versus squared transport-map error on the two-dimensional banana prior for one- and six-step iMF generation. Points represent training checkpoints, connected in training order; stars mark the final checkpoint. Error bars show one Monte Carlo standard error, conditional on the fitted model and numerical reference, and do not represent variability across training runs. B.2 Synthetic experiments Prior training. For each scenario, all models are trained on the same fixed dataset of 20,000 prior samples using four-hidden-layer SiLU networks of width 512, with 16,000 optimizer updates at learning rate 10−310^-3. The iMF model uses raw time inputs (r,t)(r,t) and separate average-velocity and auxiliary-velocity heads; FM and diffusion use four Fourier time harmonics. Both iMF and FM use independent data-noise coupling. For iMF, we sample r,tr,t independently from Unif[0,1]Unif[0,1], and set the auxiliary-loss weight to λ=1λ=1. We use Adam for iMF with weight decay 10−510^-5 and AdamW for FM and diffusion with zero weight decay, both with (β1,β2)=(0.9,0.999)( _1, _2)=(0.9,0.999). Gradient clipping at norm 5 is applied only to FM and diffusion. Figure 4 shows that all three models capture the main geometric structure of the banana and sine priors. (a) Banana prior. (b) Sine prior. Figure 4: Learning the priors by iMF, FM, and diffusion model. Posterior sampling. The two iMF samplers share a six-step transport with uniformly spaced times, while FM+SPT+pCN uses 100 RK4 steps (400 network evaluations per transport). All SPT methods use ten parallel ladders with 24 inverse temperatures βℓ=(ℓ/23)2 _ =( /23)^2, ℓ=0,…,23 =0,…,23, initialized from (0,I32)N(0,I_32). The pCN variants update every replica with pCN; the hybrid replaces the coldest update with split HMC using six integration steps and initial step size 0.050.05. During 500 adaptation sweeps, proposal parameters are shared across ladders and updated separately at each temperature using the observed acceptance fraction across the ten ladders. We apply additive updates to the logit pCN noise scale and log HMC step size, proportional to the difference between observed and target acceptance, with respective gains 0.10/1+s/300.10/ 1+s/30 and 0.055/1+s/300.055/ 1+s/30 at sweep s=0,…,499s=0,…,499. We target acceptance rates of 0.300.30 for pCN and 0.750.75 for split HMC. These are practical choices consistent with the pCN tuning range 0.30.3–0.50.5 used by Chen et al. (2016) and the general (leapfrog) HMC range 0.60.6-0.90.9 discussed by Betancourt et al. (2014). The hottest pCN noise scale remains fixed at 0.9990.999. Each sweep includes odd and even adjacent-swap passes, with their order alternating across sweeps. After adaptation, all proposal parameters are frozen, 500 additional burn-in sweeps are discarded, and 600 consecutive cold states per ladder are retained without thinning, yielding 6,000 samples per method and seed. DPS generates the same number of samples using 1,000 DDPM steps with a linear noise schedule from 10−410^-4 to 2×10−22× 10^-2. Scaling DPS guidance by σy−2 _y^-2 caused numerical instability in our experiments. Following the residual-dependent step-size scheme of Chung et al. (2023, Appendix D.1) and its official implementation, we use guidance −0.3∇xt‖G(x^0(xt))−y‖2-0.3 _x_t\|G( x_0(x_t))-y\|_2, without explicit noise-variance scaling. Sampling is repeated with seeds 0,1,2,30,1,2,3; full implementation settings are provided with the released code. Figures 5 and 6 compare sample histograms with numerical posterior references. Both iMF- and FM-based SPT recover the main posterior well, while DPS exhibits systematic deviations from the reference distributions. Evaluation metrics. Let d be the dimension, C the number of ladders, and N the number of retained samples per ladder. Writing ρj,c(ℓ) _j,c( ) for the lag-ℓ autocorrelation of coordinate j in ladder c, the mean coordinate ESS is Coord.ESS=1d∑j=1d∑c=1CN1+2∑ℓ=1∞ρj,c(ℓ),Coord.ESS/s=Coord.ESSTsampling.Coord.\ ESS= 1d _j=1^d _c=1^C N1+2 _ =1^∞ _j,c( ), .\ ESS/s= Coord.\ ESST_sampling. For distributions μ,νμ,ν on ℝdR^d, the sliced Wasserstein distance is SW2(μ,ν)=(∫d−1W22((Pω)#μ,(Pω)#ν)σ(ω))1/2,Pω(x)=ω⊤x,SW_2(μ,ν)= ( _S^d-1W_2^2 ((P_ω)_\#μ,(P_ω)_\#ν )\,dσ(ω) )^1/2, P_ω(x)=ω x, where σ is the uniform probability measure on the unit sphere. For J two-dimensional blocks and a common partition Bj,bb=1B\B_j,b\_b=1^B of each block’s first-coordinate marginal, define pj,b=μx:xj,1∈Bj,bp_j,b=μ\x:x_j,1∈ B_j,b\ and qj,b=νx:xj,1∈Bj,bq_j,b=ν\x:x_j,1∈ B_j,b\. Then BinTV(μ,ν)=1J∑j=1J12∑b=1B|pj,b−qj,b|.Bin\ TV(μ,ν)= 1J _j=1^J 12 _b=1^B|p_j,b-q_j,b|. Figure 5: Posterior sample histograms on the banana prior. Figure 6: Posterior sample histograms on the sine prior. B.3 CLIP-guided ImageNet experiments B.3.1 Experiment details We use the pretrained iMF-XL/244 4 https://github.com/Lyy-iiis/imeanflow/tree/torch with one transport step, ω=8ω=8, and guidance interval [0.4,0.65][0.4,0.65], and frozen SiT-XL/255 5 https://huggingface.co/sairights/sit-xl-2-256x256-sde-cfg/resolve/5a15dcc/pretrained_models/SiT-XL-2-256x256.safetensors with 125 Heun steps (250 velocity evaluations) and classifier-free guidance scale 44. For each class-prompt pair, we use eight independent ladders of eight replicas with linear inverse temperatures βk=k/7 _k=k/7, initialized from (0,I)N(0,I). The iMF methods run 250 adaptation sweeps, 250 burn-in sweeps, and 1,000 retained sweeps per ladder, yielding 8,000 retained cold states. At each adaptation sweep, we average the local proposal acceptance rates over all eight ladders. We then make progressively smaller updates to the pCN proposal scale to target 0.300.30 acceptance and, for iMF+SPT+hybrid, to the split-HMC step size at the coldest replica to target 0.750.75 acceptance. The resulting proposal parameters are fixed before burn-in. iMF+SPT+hybrid uses six-step split HMC at the coldest replica and pCN at all remaining replicas. SiT+SPT+pCN is run after the two iMF methods under a prompt-specific budget: its sampling budget equals the measured iMF+SPT+hybrid sampling time, with an additional one-hour allowance for loading, preflight, scoring, and export. Its available sampling time is allocated as 1/61/6, 1/61/6, and 2/32/3 to adaptation, burn-in, and retained sampling; only complete sweeps are performed. B.3.2 Full results Table 3 reports the complete results, and Figures 8-18 display samples from the first two ladders for each method. Table 3: CLIP reward over retained samples and sampling cost. Method Mean CLIP reward ↑ Sampling time (h) Mean sweep time (s) Golden retriever: “A photograph of a golden retriever running through fallen autumn leaves.” Best-of-K 0.2864 2.49 – SiT+SPT+pCN 0.3122 10.51 370.89 iMF+SPT+pCN 0.3809 2.40 5.75 iMF+SPT+hybrid 0.3824 10.95 26.28 Golden retriever: “A beautiful professional photograph of a golden retriever in a sunlit meadow, sharp focus, natural colors.” Best-of-K 0.3073 2.48 – SiT+SPT+pCN 0.3368 10.79 373.28 iMF+SPT+pCN 0.3741 2.43 5.83 iMF+SPT+hybrid 0.3754 11.26 27.02 Golden retriever: “A colorful painting of a golden retriever standing in falling snow.” Best-of-K 0.2885 2.45 – SiT+SPT+pCN 0.3355 10.96 372.17 iMF+SPT+pCN 0.3945 2.50 6.01 iMF+SPT+hybrid 0.4011 11.50 27.58 Golden retriever: “A photograph of a golden retriever swimming in a clear blue lake.” Best-of-K 0.3064 2.45 – SiT+SPT+pCN 0.3675 10.69 369.97 iMF+SPT+pCN 0.3847 2.45 5.88 iMF+SPT+hybrid 0.3845 11.16 26.69 Bald eagle: “A photograph of a bald eagle perched on a branch surrounded by red and golden autumn leaves.” Best-of-K 0.3086 2.50 – SiT+SPT+pCN 0.3156 10.73 371.46 iMF+SPT+pCN 0.3743 2.46 5.90 iMF+SPT+hybrid 0.3781 11.20 26.87 Bald eagle: “A photograph of a bald eagle catching a fish from a lake, with splashing water.” Best-of-K 0.3247 2.59 – SiT+SPT+pCN 0.3407 10.78 376.60 iMF+SPT+pCN 0.3625 2.46 5.90 iMF+SPT+hybrid 0.3596 11.17 26.81 Bald eagle: “A photograph of a bald eagle soaring above snow-capped mountains, wings spread wide.” Best-of-K 0.3278 2.49 – SiT+SPT+pCN 0.3481 10.70 370.27 iMF+SPT+pCN 0.3936 2.44 5.85 iMF+SPT+hybrid 0.3914 11.15 26.76 Bald eagle: “A colorful painting of a bald eagle perched on a pine branch in falling snow.” Best-of-K 0.2984 2.58 – SiT+SPT+pCN 0.3139 10.74 371.58 iMF+SPT+pCN 0.3727 2.41 5.77 iMF+SPT+hybrid 0.3713 11.16 26.78 Mink: “A photograph of a mink beside a clear forest stream, surrounded by green moss.” Best-of-K 0.3629 2.75 – SiT+SPT+pCN 0.3946 10.72 371.04 iMF+SPT+pCN 0.4190 2.45 5.89 iMF+SPT+hybrid 0.4112 11.12 26.70 Mink: “A colorful painting of a mink standing in falling snow.” Best-of-K 0.3160 2.49 – SiT+SPT+pCN 0.3487 10.47 369.64 iMF+SPT+pCN 0.4186 2.36 5.66 iMF+SPT+hybrid 0.4162 10.99 26.37 Castle: “A colorful painting of a medieval castle covered in snow.” Best-of-K 0.2489 2.48 – SiT+SPT+pCN 0.2797 10.53 371.82 iMF+SPT+pCN 0.3618 2.42 5.79 iMF+SPT+hybrid 0.3569 10.92 26.21 Castle: “A photograph of a medieval castle on a hill at sunrise, surrounded by green countryside.” Best-of-K 0.3023 2.46 – SiT+SPT+pCN 0.3476 10.96 372.13 iMF+SPT+pCN 0.3820 2.52 6.05 iMF+SPT+hybrid 0.3849 11.40 27.35 Aircraft carrier: “A photograph of an aircraft carrier sailing across the ocean at sunset.” Best-of-K 0.2879 2.49 – SiT+SPT+pCN 0.3286 10.38 373.72 iMF+SPT+pCN 0.3667 2.42 5.80 iMF+SPT+hybrid 0.3659 10.90 26.16 Aircraft carrier: “A watercolor painting of an aircraft carrier on a calm blue sea.” Best-of-K 0.2956 2.52 – SiT+SPT+pCN 0.3055 10.63 367.95 iMF+SPT+pCN 0.3736 2.49 5.97 iMF+SPT+hybrid 0.3774 11.07 26.55 (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 7: Class: golden retriever, prompt: “A beautiful professional photograph of a golden retriever in a sunlit meadow, sharp focus, natural colors.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 8: Class: bald eagle, prompt: “A photograph of a bald eagle perched on a branch surrounded by red and golden autumn leaves.” (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 9: Class: golden retriever, prompt: “A colorful painting of a golden retriever standing in falling snow.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 10: Class: bald eagle, prompt: “A colorful painting of a bald eagle perched on a pine branch in falling snow.” (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 11: Class: golden retriever, prompt: “A photograph of a golden retriever swimming in a clear blue lake.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 12: Class: bald eagle, prompt: “A photograph of a bald eagle catching a fish from a lake, with splashing water.” (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 13: Class: mink, prompt: “A photograph of a mink beside a clear forest stream, surrounded by green moss.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 14: Class: mink, prompt: “A colorful painting of a mink standing in falling snow.” (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 15: Class: castle, prompt: “A photograph of a medieval castle on a hill at sunrise, surrounded by green countryside.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 16: Class: castle, prompt: “A colorful painting of a medieval castle covered in snow.” (a) best-of-K. (b) SiT + SPT + pCN. (c) iMF + SPT + pCN. (d) iMF + SPT + hybrid. Figure 17: Class: aircraft carrier, prompt: “A photograph of an aircraft carrier sailing across the ocean at sunset.” (e) best-of-K. (f) SiT + SPT + pCN. (g) iMF + SPT + pCN. (h) iMF + SPT + hybrid. Figure 18: Class: aircraft carrier, prompt: “A watercolor painting of an aircraft carrier on a calm blue sea.”