Tail Annealing for Heavy-Tailed Flow Matching
Jean Pachebat 1
arXiv:2605.20068v1 [stat.ML] 19 May 2026
Abstract
both in theory and in practice. Jaini et al. (2020) proved that Lipschitz transformations preserve tail type: light-tailed inputs yield light-tailed outputs. Since most practical architectures (MLPs in GANs, coupling layers in normalizing flows, and diffusion denoisers) use Lipschitz-continuous functions for training stability, they cannot generate heavy tails from Gaussian noise without explicit tail-modifying components.
Standard generative models struggle with heavytailed data: Lipschitz architectures cannot produce power-law tails from Gaussian noise, and interpolating between heavy-tailed data and Gaussians is ill-posed. We propose a simple fix: apply the soft-log transform ϕ(x) = sign(x) · log(1 + |x|) coordinate-wise to data before training, then exponentiate samples after generation. A Hill diagnostic decides per-coordinate whether to transform, leaving light-tailed margins untouched at no added complexity. This compresses heavy tails into a range where standard flow matching succeeds, without heavy-tailed base distributions or architectural modifications. We provide theoretical intuition for why this works: the log-transform maps Pareto tails to exponentials, and the induced dynamics implement a form of tail annealing via power transformations. On a 144-configuration multivariate benchmark (3 copulas, d up to 100, 4 tail indices), Log-FM dominates specialized baselines on W1 , CVaR99 , and extreme-quantile metrics, and is the only method with zero severe divergences across 2,880 runs.
Prior approaches. Existing methods address heavy tails through two strategies. Jaini et al. (2020) replace Gaussian noise with heavy-tailed base distributions (Student-t), but this introduces training instability. Hickling & Prangle (2025) take a different approach: generate samples from a standard flow, then apply a tail-modifying transformation to the output. Their method requires estimating the tail parameter λ = 1/ν via the Hill estimator. Our approach. We propose a simple fix: transform data via ϕ(x) = sign(x) · log(1 + |x|) before training, run standard flow matching in log-space, then apply ϕ−1 to generated samples. Unlike Hickling & Prangle (2025) who transform outputs with estimated tail parameters, our input transform is parameter-free: the transformation ϕ works regardless of the true tail index. We provide theoretical intuition for why this works. The logtransform compresses heavy tails: Pareto becomes approximately exponential. In log-space, both data and Gaussian noise have light tails, so standard interpolation is well-posed. The induced dynamics in original space implement a form of tail annealing via power transformations X0αt .
1. Introduction Heavy-tailed distributions are ubiquitous in real-world phenomena: financial returns, insurance claims, earthquake magnitudes, network degrees, and climate extremes. In these domains, rare events dominate outcomes: asset returns exhibit empirically documented fat tails that Gaussian models systematically underestimate (Cont, 2001), and similar departures from Gaussianity are well-known in insurance (Embrechts et al., 2013) and natural-hazard data. Generating realistic synthetic data for stress testing, simulation, or data augmentation requires generative models that capture tail behavior.
Contributions. • We show that a parameter-free log-transform makes standard flow matching work for heavy-tailed data, without requiring tail estimation or architectural changes (§4). • We analyze the tail annealing mechanism: the logtransform maps heavy tails to exponential tails, and the induced dynamics implement power transformations X0αt that continuously adjust the tail index (§3).
Yet standard generative models fail on heavy-tailed data, 1
CMAP, École Polytechnique, Institut Polytechnique de Paris, France. Correspondence to: <[email protected]>.
• Our 144-configuration multivariate benchmark (3 copulas, 4 dimensions up to d = 100, 4 tail indices, 20
Preprint. May 20, 2026.
1
Tail Annealing for Heavy-Tailed Flow Matching
Diffusion and flow matching. Denoising diffusion models (Ho et al., 2020; Song et al., 2021b) and flow matching (Lipman et al., 2023) construct a path from data to noise and learn to reverse it. Given data X0 and noise X1 ∼ N (0, I), the interpolant is:
replications, two margin types) shows Log-FM dominates baselines on tail metrics and is the only method with zero severe divergences (§5). Real-data validation on Fama–French is in Appendix E.2.
2. Background
Xt = αt X0 + βt X1 ,
We review heavy-tailed distributions and the fundamental barriers they pose for generative modeling.
with schedules satisfying (α0 , β0 ) = (1, 0) and (α1 , β1 ) = (0, 1). A neural network learns the velocity field vθ (xt , t) or score ∇ log pt (xt ), enabling generation by integrating from noise to data.
2.1. Heavy-Tailed Distributions A distribution is heavy-tailed if its tails decay slower than any exponential, formally E[eλX ] = ∞ for all λ > 0 (Nair et al., 2022). Heavy-tailed distributions arise naturally in finance (Cont, 2001), insurance, climate science, and network traffic, where rare events have outsized impact.
2.3. The Lipschitz Barrier Jaini et al. (2020) proved that Lipschitz maps preserve tail type: Theorem 2.1 (Jaini et al. (2020)). If Z is light-tailed and T is Lipschitz, then T (Z) is light-tailed.
The canonical example is the Pareto distribution with survival function P(X > t) = t−1/γ for t ≥ 1, where γ > 0 is the shape parameter (the GPD shape, equivalently the extreme-value index). Larger γ means heavier tails: γ < 1 for finite mean, γ < 1/2 for finite variance. The Student-t distribution with ν degrees of freedom has Pareto-like tails with effective γ = 1/ν. Lognormal has γ = 0 (Gumbel max-domain of attraction): heavy in the MGF sense, but with no polynomial decay rate.
Neural networks composed of linear layers and standard activations (ReLU, tanh, sigmoid) are Lipschitz. Thus, standard normalizing flows cannot map Gaussian noise to heavytailed outputs. 2.4. Prior Approaches Heavy-tailed base distributions. Tail-Adaptive Flows (TAF; Jaini et al., 2020) use Student-t base distributions. Extensions include marginal TAF (Laszkiewicz et al., 2022) and generalized TAF. Heavy-tailed diffusion (Pandey et al., 2024) replaces Gaussian noise with Student-t. These methods address the Lipschitz barrier but can suffer training instability.
The log-transform. A key observation: the logarithm maps heavy tails to lighter tails. If X ∼ Pareto(γ), then log X ∼ Exp(1/γ), an exponential distribution with rate 1/γ. This reflects that Pareto belongs to an exponential family with sufficient statistic log x. The log-transform thus provides a natural bridge between power-law and exponentialfamily distributions; we develop the formal theory in Section 3.
Tail-modifying transforms. Tail Transform Flows (TTF; Hickling & Prangle, 2025) add explicit layers that map Gaussian tails to power-law tails using erfc-based transformations. This requires estimating tail parameters per dataset via the Hill estimator.
Regular variation. Heavy-tailed distributions are characterized by regular variation: P(X > tx)/P(X > t) → x−α as t → ∞, with tail index α > 0. This captures the essential property of polynomial tail decay without requiring the exact Pareto form. For X ∼ Pareto(γ), the corresponding tail index is α = 1/γ; Student-t, Burr, and log-gamma are also regularly varying. The log-transform maps any regularly varying distribution to one with exponential-type tails; see Proposition 3.8.
In the next section, we present an alternative: apply the logtransform to data before training, work entirely in log-space, then exponentiate samples after generation. This parameterfree approach requires no tail estimation and we show that it circumvents the Lipschitz barrier.
2.2. Generative Models
3. Tail Annealing in Log-Space
Normalizing flows. Normalizing flows (Rezende & Mohamed, 2015; Dinh et al., 2017) learn an invertible transformation T : Rd → Rd mapping a base distribution (typically Gaussian) to the target. The change-of-variables formula gives exact likelihoods. However, for stability, flow architectures use Lipschitz-continuous components.
This section develops the theoretical foundation for logspace flow matching. The log-transform compresses heavy tails into a light-tailed regime where standard methods apply, and the induced dynamics in original space implement tail annealing: power transformations that continuously adjust the tail index from heavy to light. 2
Tail Annealing for Heavy-Tailed Flow Matching
3.2. Induced Process in Original Space
The principle of annealing complex structure through a sequence of easier problems is well-established in machine learning. Noise annealing in score-based models (Song & Ermon, 2019) gradually reduces noise levels to enable learning across scales. Curriculum learning (Bengio et al., 2009) trains on easy examples before hard ones. Our tail annealing follows the same principle: rather than directly modeling heavy-tailed data (hard), we work in log-space where tails are light (easy), and the power transformation X0αt provides a principled path through intermediate tail weights. Section 4 presents the complete algorithm.
We apply flow matching to the transformed data X̃0 = ϕ(X0 ), interpolating with Gaussian noise X̃1 ∼ N (0, Id ) via X̃t = αt X̃0 + βt X̃1 (see Section 4 for the complete framework). Here we analyze the induced process in original space. Applying ϕ−1 to the interpolant yields, for large positive values where ϕ ≈ log: Xt := ϕ−1 (X̃t ) ≈ X0αt · eβt X̃1 .
(1)
This factorizes into a power-transformed data term X0αt and a log-normal noise term eβt X̃1 . Although eβt X̃1 is itself heavy-tailed in the MGF sense, it has extreme-value index 2 γ = 0 (Gumbel MDA): P(Yt > t) ∼ e−c(log t) decays faster than any power law and all moments are finite. By Breiman’s lemma it does not contribute a polynomial tail of its own to the product:
3.1. The Soft-Log Transform We work with the soft-log transform ϕ : R → R defined componentwise by: ϕ(x) = sign(x) · log(1 + |x|); ϕ−1 (y) = sign(y) · (e|y| − 1).
Proposition 3.1 (Product Preserves Regular Variation; Breiman). Let X be regularly varying with index −α (α > 0), and let Y > 0 be independent with E[Y p ] < ∞ for all p > 0 (e.g., any variable in the Gumbel MDA, such as log-normal). Then XY is regularly varying with index −α.
Unlike the standard logarithm, ϕ is smooth at the origin (ϕ′ (0) = 1) and defined on all of R. For large |x|, we have ϕ(x) ≈ sign(x) log |x|, so tail behavior matches the logarithm. Relation to Box-Cox. Transforming heavy-tailed data through a log transform is a known trick in the statistical litterature. The soft-log is inspired by, but distinct from, the Box-Cox family (Box & Cox, 1964). Recall Box-Cox: ϕλ (x) = (xλ − 1)/λ for λ ̸= 0, and ϕ0 (x) = log x for x > 0. Our transform differs in two ways: it extends to all of R via the signed construction, and uses log(1+|x|) rather than log |x| for smoothness at zero. Asymptotically, both behave as log |x| for large |x|. Methods like TTF (Hickling & Prangle, 2025) use power transforms with λ = 1/ν estimated via the Hill estimator. In our approach, numerical computation make no use of numerical values of the tail estimate, which are notoriously unstable: the same transform applies regardless of the true tail index.
This classical result (see Resnick (1987)) ensures that the induced process Xt inherits its polynomial-tail behavior from X0αt , not from eβt X̃1 (whose extreme-value index is 0). The power transformation is therefore what drives tail annealing.
The soft-log maps heavy tails to light tails. If X ∼ Pareto(γ) with P(X > t) = t−1/γ (so γ > 0 is the shape parameter and 1/γ is the tail index), then for large y:
1 ∇x̃ log pX̃ (x̃) ≈ − . γ
3.3. Score Functions in Log-Space The score function ∇x̃ log pt (x̃) in log-space has fundamentally different behavior than in original space. Proposition 3.2 (Log-Space Score of Pareto). Let X ∼ Pareto(γ) (so F̄X (t) = t−1/γ ). Then X̃ = ϕ(X) has approximately exponential tails: for large x̃,
The score is approximately constant in the tails.
P(ϕ(X) > y) = P(X > ey − 1) ≈ e−y/γ ,
Proof. For large x, ϕ(x) ≈ log x, so X̃ ≈ log X ∼ Exp(1/γ). The exponential density is p(x̃) ∝ e−x̃/γ , giving ∇x̃ log p(x̃) = −1/γ.
i.e., X̃ = ϕ(X) has exponential-type tails. More generally, for any distribution in the Fréchet domain of attraction, ϕ yields a distribution with exponential or lighter tails.
Pareto score. By the Pareto score we mean the score of the Pareto density pX (x) = (1/γ) x−1/γ−1 in original space:
Since ϕ is a diffeomorphism, densities transform via the standard change of variables: pX̃ (x̃) = pX (ϕ−1 (x̃)) · |(ϕ−1 )′ (x̃)| = pX (ϕ−1 (x̃)) · e|x̃| .
sX (x) = ∇x log pX (x) = − 3
1/γ + 1 , x
Tail Annealing for Heavy-Tailed Flow Matching
which diverges as x → 0+ and decays only as 1/x in the tails. The log-space score in Proposition 3.2 is by contrast bounded by 1/γ everywhere, so the velocity field that flow matching has to regress is Lipschitz in x̃. This score boundedness is what enables the standard MLP velocity network to learn a target with regularly varying tails in the original space.
Corollary 3.6 (Tail Lightening Along Forward Process). Along the forward process with αt decreasing from 1 to 0: • At t = 0: α0 = 1; X0α0 = X0 ∼ Pareto(γ) (original heavy tails) • At intermediate t: X0αt ∼ Pareto(γαt ) (lighter tails as αt → 0)
Proposition 3.3 (Score of Noised Log-Data). For the interpolant X̃t = αt X̃0 + βt X̃1 with X̃0 = ϕ(X0 ) having exponential-type tails and X̃1 ∼ N (0, Id ), the marginal score satisfies: ∇x̃t log pt (x̃t ) = −
• As αt → 0: X0αt → 1 almost surely (degenerate) The full process Xt = X0αt ·eβt X̃1 transitions from regularlyvarying Pareto (shape γ, tail index 1/γ) to log-normal (extreme-value index 0, no polynomial tail).
ˆ1 (x̃t , t) x̃ , βt
This is what we frame as the tail annealing mechanism: the power transformation X0αt preserves the Pareto family while continuously adjusting the shape parameter from γ (heavy) toward 0 (light) as αt → 0. We state Lemma 3.4 for α ∈ (0, ∞) rather than the forward-process range αt ∈ (0, 1], because the same identity appears in Proposition 3.9 for the regularly varying extension.
ˆ1 (x̃t , t) := E[X̃1 | X̃t = x̃t ] is the conditional where x̃ expectation of the noise. This is the standard score-denoiser relationship. The key difference from original-space diffusion is that the transformed data distribution pX̃0 is light-tailed, so: 1. The score ∇ log pt is well-behaved for all t
3.5. Extension to Regularly Varying Distributions
ˆ1 (x̃t , t) has bounded outputs 2. The denoiser x̃
The Pareto-exponential correspondence extends to the broader class of regularly varying distributions.
3. Standard neural network architectures suffice
Definition 3.7 (Regularly Varying). A cumulative distribution function F on R is regularly varying with index −α (written F ∈ RV−α ) if its survival function F̄ (x) = 1 − F (x) satisfies:
3.4. Tail Behavior Under Power Transformation The induced process (1) involves the power-transformed data X0αt . We characterize its distribution. Note that this analysis describes the decoded interpretation: in practice, all computation occurs in log-space where distributions are light-tailed. The heavy-tailed structure below is what we would observe if we applied ϕ−1 to the log-space interpolant at each t.
lim
t→∞
F̄ (tx) = x−α F̄ (t)
for all x > 0. Here α > 0 is the (polynomial) tail index, in the standard EVT convention; for X ∼ Pareto(γ), α = 1/γ. Regular variation captures the essential property of polynomial tail decay without requiring the exact Pareto form:
Lemma 3.4 (Power Transformation Preserves Pareto Family). Let X0 ∼ Pareto(γ) with density pX0 (x) = (1/γ) x−1/γ−1 for x ≥ 1, so γ > 0 is the shape parameter (larger γ means heavier tails). For any α ∈ (0, ∞), X0α ∼ Pareto(γα): the shape parameter scales linearly with the exponent.
Distribution Pareto(γ) Student-t(ν) Burr(c, k) Log-gamma
Proof. The survival function of Y = X0α is P(Y > y) = P(X0 > y 1/α ) = y −1/(γα) for y ≥ 1, which is the survival function of Pareto(γα).
Tail behavior −1/γ
Index
F̄ (x) = x F̄ (x) ∼ c · x−ν F̄ (x) ∼ c′ · x−ck F̄ (x) ∼ c′′ · x−α (log x)β
−1/γ −ν −ck −α
Proposition 3.8 (Log-Transform of Regularly Varying; informal). If X is regularly varying with index −α, then X̃ = ϕ(X) has exponential-type right tail with rate α:
Remark 3.5 (Parametrization). We use γ for the GPD shape parameter / extreme-value index, so P(X > t) = t−1/γ . The corresponding tail index (Hill estimator output, the standard EVT parameter) is 1/γ. Under Lemma 3.4 the power transform X0α has shape parameter γα and tail index 1/(γα).
− log P(X̃ > z) = α z + o(z)
as z → ∞.
A precise two-sided Potter-bound version is stated and proved in Appendix C (Proposition C.1). The intuition is direct: writing the regularly varying tail as F̄X (x) = 4
Tail Annealing for Heavy-Tailed Flow Matching
x−α L(x) with L slowly varying, and using ϕ−1 (z) = ez − 1 ∼ ez , gives P(X̃ > z) = e−αz L(ez ), and the slowly varying factor only contributes a subexponential correction log L(ez ) = o(z).
Heterogeneous margins. When the coordinates of X have different tail behaviour (heavy-tailed in some, lighttailed in others), applying ϕ to every coordinate slightly distorts the light-tailed margins. The adaptive variant of Section 4.2 addresses this by Hill-gating the transform coordinate-wise; a continuous generalization of the soft-log via a scale parameter s2 , of which the adaptive variant is the binary s2 ∈ {0, 1} instance, is developed in Appendix D.
This proposition is the conceptual reason flow matching in log-space is well-posed for the full class of regularly varying distributions, not just exact Pareto. Up to slowly-varying corrections, X̃ has the same exponential rate e−αz as log X, which equals e−z/γ when X ∼ Pareto(γ). The score in logspace is therefore bounded in the tails and standard diffusion methods apply.
3.7. Circumventing the Lipschitz Barrier The Lipschitz barrier of Jaini et al. (2020) states that Lipschitz maps preserve tail type: a normalizing flow with Lipschitz layers cannot map Gaussian noise to heavy-tailed outputs. Our construction circumvents this by placing the non-Lipschitz transformation (ϕ−1 ) outside the learned dynamics. The flow operates entirely in log-space where both endpoints (exponential-type transformed data and Gaussian noise) have light tails. Heavy tails emerge only upon applying ϕ−1 at the final sampling step.
Beyond regular variation. The transform is also beneficial for subexponential distributions outside the Fréchet domain. For Weibull tails F̄ (x) ∼ exp(−xβ ) with β < 1, the log-transform yields a doubly-exponential tail P(X̃ > z) ∼ exp(−(ez − 1)β ), which is well within the regime where standard flow matching has no difficulty. For lognormal data, ϕ(X) ≈ log X is exactly Gaussian, so the log-space target is the ideal case for FM. Light-tailed margins are addressed by the adaptive variant of Section 4.2.
4. Algorithm
Proposition 3.9 (Power Transformation of Regularly Varying). If X ∈ RV−α , then X β for β ∈ (0, ∞) is regularly varying with index −α/β.
We present the complete algorithm for log-space generative modeling. The method requires the three following modifications to standard flow matching: (1) transform data via ϕ before training (with a Hill-based diagnostic deciding which coordinates to transform), (2) train a velocity network on the transformed variables, and (3) apply ϕ−1 to the corresponding coordinates of generated samples. The Hill diagnostic distinguishes Log-FM from a standard coordinate-wise logtransform: it lets the method handle heterogeneous-margin data with a single Hill estimate per coordinate at fit time.
Proof. For Y = X β , we have P(Y > ty)/P(Y > t) = P(X > (ty)1/β )/P(X > t1/β ) → y −α/β as t → ∞, by regular variation of X. This is the regularly varying analogue of Lemma 3.4: power transformations remap the tail index, so X αt along the forward process continuously interpolates between the original tail (αt = 1) and arbitrarily light tails (αt ∈ [0, 1)]).
4.1. Flow Matching Framework
This confirms that the tail-annealing mechanism of Lemma 3.4 extends beyond exact Pareto to the full class of regularly varying distributions.
We use the flow matching framework (Lipman et al., 2023), which can also be viewed through the lens of stochastic interpolants (Albergo & Vanden-Eijnden, 2023).
3.6. Multivariate Extension Forward process. Given log-transformed data X̃0 = ϕ(X0 ) and noise X̃1 ∼ N (0, Id ), define the interpolant:
In the multivariate setting, X ∈ Rd , we apply ϕ coordinatewise: ϕ(x) = (ϕ(x1 ), . . . , ϕ(xd )),
X̃t = αt X̃0 + βt X̃1 ,
X̃ = ϕ(X).
where (αt , βt ) are differentiable schedules satisfying boundary conditions (α0 , β0 ) = (1, 0) and (α1 , β1 ) = (0, 1). At t = 0, X̃0 is the transformed data; at t = 1, X̃1 is pure Gaussian noise. The conditional distribution is Gaussian: qt|0 (x̃t | x̃0 ) = N (x̃t ; αt x̃0 , βt2 Id ).
Each coordinate of X̃ inherits the exponential-rate tail of Proposition 3.8, so the marginals of X̃ all live in the Gumbel domain. Because ϕ is a diffeomorphism with diagonal Jacobian Jϕ (x) = diag((1 + |xi |)−1 ), the dependence structure (copula) of X is preserved exactly by X̃: the log-transform is a marginal operation. Interpolating X̃t = αt X̃0 + βt X̃1 with X̃1 ∼ N (0, Id ) is therefore well-posed even when X has strong extremal dependence: the velocity network sees only the light-tailed transformed variates and never has to extrapolate to the polynomial regime.
Velocity field. Differentiating the interpolant gives the per-trajectory time derivative X̃˙ t = α̇t X̃0 + β̇t X̃1 , 5
Tail Annealing for Heavy-Tailed Flow Matching
Algorithm 1 Log-FM: Training
where α̇t = dαt /dt and β̇t = dβt /dt. The corresponding marginal velocity field, vt (x̃t ) := E[α̇t X̃0 + β̇t X̃1 | X̃t = x̃t ], is what a neural network vθ regresses against the pertrajectory target:
d Input: Dataset {xi }N i=1 ⊂ R , velocity network vθ , schedule (αt , βt ), iterations T , Hill threshold αmax = 4 Fit transform: for each coordinate j, compute Hill esti(j) mate α̂j on {xi }i ; set mask mj via (2) Preprocess: x̃i ← Φ(xi ) for all i {coordinate-wise; identity where mj = 0} for iter = 1 to T do (j) Sample batch {x̃0 } from {x̃i } (j) Sample x̃1 ∼ N (0, Id ) Sample t ∼ Uniform[0, 1] (j) (j) (j) Form interpolant x̃t ← αt x̃0 + βt x̃1 (j) (j) Compute target u(j) ← α̇t x̃0 + β̇t x̃1 P (j) 1 (j) 2 Loss L ← |batch| j ∥vθ (x̃t , t) − u ∥ Update θ via gradient descent on L end for Output: Trained vθ , mask (mj )dj=1
L(θ) = Et∼U [0,1] EX̃0 ,X̃1 2 vθ (X̃t , t) − (α̇t X̃0 + β̇t X̃1 ) . Connection to score and denoiser. The velocity field reˆ0 (x̃t , t) := E[X̃0 | lates to the score ∇ log pt and denoiser x̃ X̃t = x̃t ] via: ˆ0 (x̃t , t) + β̇t x̃ ˆ1 (x̃t , t) vt (x̃t ) = α̇t x̃ ˆ1 (x̃t , t) x̃ ∇x̃t log pt (x̃t ) = − , βt ˆ1 = (x̃t − αt x̃ ˆ0 )/βt . See Appendix B for derivawhere x̃ tions.
with inverse Φ−1 (x̃)j = ϕ−1 (x̃j ) when mj = 1 and x̃j otherwise. The threshold αmax = 4 is not sensitive: on heterogeneous-margin data with Pareto coordinates (α̂ ≈ 1.5–2.5) and Gaussian coordinates (α̂ ≈ 6), any threshold in [3, 5] yields the same mask. Setting mj ≡ 1 recovers the uniform-ϕ variant analyzed in Section 3; the adaptive rule is the discrete instance of a continuous parametrized soft-log (j) family φs2 with s2 ∈ {0, 1} chosen by the Hill diagnostic, developed in Appendix D. We use Φ (with the Hill-based mask) as our default; we also report the uniform variant as an ablation in Section 5.3.
Schedule choices. The interpolation schedule (αt , βt ) controls the path from data to noise. Common choices include the linear schedule (optimal transport paths), variance-preserving (VP) schedules, and quadratic schedules that anneal tails more aggressively; see Table 6 in Appendix B. We use the linear schedule throughout; preliminary tests showed minimal sensitivity to schedule choice, consistent with prior flow matching work (Lipman et al., 2023). Connection to tail annealing. Recall from Section 3 that the induced process in original space is approximately Xt ≈ X0αt · eβt X̃1 . The schedule αt thus controls the rate of tail annealing: faster decay of αt (e.g., quadratic) anneals tails more aggressively early in the process, while slower decay (e.g., VP polynomial) preserves heavier tails longer.
4.3. Training Algorithm 1 summarizes the training procedure. The only departure from standard flow matching is the initial logtransform of data (with the Hill diagnostic of Section 4.2).
4.2. Adaptive Coordinate Selection
4.4. Sampling
Applying the soft-log to a light-tailed coordinate compresses its bulk for no benefit. We therefore select per coordinate whether to transform, using a Hill estimate on the training marginals. Let α̂j be the Hill estimator on the upper order (j) statistics of {xi }N i=1 (an estimate of the polynomial tail index, with α̂j = 1/γ̂j under the Pareto-shape parametrization). Define the mask
Generation reverses the flow: integrate the learned velocity field from noise (t = 1) to data (t = 0), then apply the inverse transform. The velocity field satisfies the ODE:
mj = 1{α̂j ≤ αmax },
αmax = 4,
and the per-coordinate transform Φ : Rd → Rd by ( ϕ(xj ) if mj = 1, Φ(x)j = xj otherwise,
dX̃t = vθ (X̃t , t), dt which we solve backward in time. Algorithm 2 uses Euler discretization; higher-order solvers (Heun, RK4) improve sample quality with fewer steps.
(2)
4.5. Practical Considerations Output clamping (optional). The inverse transform ϕ−1 (x̃) = sign(x̃)(e|x̃| − 1) is exponential, so small errors
(3) 6
Tail Annealing for Heavy-Tailed Flow Matching
Algorithm 2 Log-FM: Sampling
4.6. Extensions
Input: Trained vθ , mask (mj ), steps K, optional clamp c (default c = ∞) Sample x̃ ∼ N (0, Id ) {initialize at t = 1} ∆t ← 1/K for k = K downto 1 do t ← k/K x̃ ← x̃ − ∆t · vθ (x̃, t) {Euler step} end for (optional) clamp coordinates with mj = 1: x̃j ← clamp(x̃j , −c, c) {only needed for α < 1} Output: x = Φ−1 (x̃) {coordinate-wise inverse from (3)}
Signed and multivariate data. The soft-log ϕ(x) = sign(x) log(1 + |x|) applies directly to signed data and acts coordinate-wise on Rd via Φ; the dependence structure is preserved exactly on transformed coordinates and learned implicitly by the velocity field. Arcsinh variant. Replacing ϕ by arcsinh(x) = log(x + √ 1 + x2 ) gives the Arcsinh-FM variant evaluated in Section 5. Both transforms have the same logarithmic asymptotic, so Propositions 3.8–3.9 apply unchanged. They differ only in regularity at the origin: arcsinh is C ∞ on R, while ϕ is C 1 but not C 2 at 0 (its second derivative jumps from +1 to −1 at the origin), which isolates whether higher-order regularity matters in practice.
deep in the tail of x̃ are exponentially amplified in x. An optional clamp x̃j ← clamp(x̃j , −c, c) before ϕ−1 provides a numerical safeguard in extremely pathological regimes only (typically α < 1: Cauchy and heavier, e.g. Hickling’s ν = 1/2 Student-t baseline). For every configuration in our main benchmark (α ≥ 1.5) the clamp is inert at c ≥ 10, and we report all main-benchmark numbers without clamping; see the ablation in Section 5.3.
5. Experiments We evaluate Log-FM against state-of-the-art heavy-tailed generative models on a controlled multivariate benchmark with known marginal tails and known copula structure. Our goals are to test whether the log-transform (i) preserves the dependence structure of the data, (ii) yields stable training across heavy-tail regimes, (iii) recovers risk-relevant tail metrics (VaR99 , CVaR99 , Q99.9 ) that are sensitive to the deep tail. A Student-t benchmark from (Hickling & Prangle, 2025) is detailed for continuity in Appendix E.1; real-data results on Fama–French are deferred to Appendix E.2.
Number of integration steps. We use K = 100 Euler steps for all experiments. Higher-order solvers (Heun, RK4) can reduce step count but Euler suffices.
5.1. Setup
Architecture. We use MLPs with sinusoidal time embeddings (each scalar t is mapped to ω [sin(2πωk t), cos(2πωk t)]K k=1 with geometric frequencies, exactly as in Ho et al. (2020)): 4 hidden layers of width 256 with SiLU activations. No architectural modifications are needed for heavy tails; the log-transform handles tail behavior.
Data. We use 144 configurations covering three copula families, three dependence strengths, four dimensions, and four tail indices: • Copulas: Gaussian copula (asymptotic independence), Gumbel copula (symmetric upper tail dependence) at Kendall’s τ ∈ {0.25, 0.5, 0.75}, and Hüsler–Reiss copula with AR(1) variogram Γij = 2(1 − ρ|i−j| ) at ρ ∈ {0.1, 0.5, 0.9}.
Likelihood evaluation. We report NLL using the continuous change-of-variables formula for continuous normalizing flows. The trace of the velocity Jacobian is estimated stochastically via Hutchinson’s trick with K = 10 Rademacher projections, and the augmented ODE is integrated by the adaptive Dormand–Prince dopri5 scheme with atol = rtol = 10−5 . The transform itself contributes P the closed-form Jacobian term j:mj =1 log(1 + |xj |), which is exact and negligible to compute.
• Margins: 70% symmetrized Pareto with tail index α ∈ {1.5, 1.75, 2.0, 2.5} (smaller = heavier), 30% standard Gaussian. • Dimensions: d ∈ {10, 20, 50, 100}. • Samples: ntrain = 10,000, nval = 5,000, ntest = 20,000.
Training. AdamW optimizer with learning rate 5 × 10−3 and weight decay 10−5 . We use full-batch gradient descent with gradient clipping at 10.0 (before inverse transform). Training runs for up to 5000 epochs with early stopping (patience 100) based on validation loss.
For each configuration we train 5 methods × 20 independent replications, giving 14,400 trained models in total. Methods. Log-FM (default, with Hill-gated transform of Section 4.2, αmax = 4); Log-FM (uniform, mj ≡ 1); 7
Tail Annealing for Heavy-Tailed Flow Matching Table 1. Median W1P across Gumbel + Gaussian copulas (lower is better, bold = best). Log-FM dominates at every tail index and dimension; the gap widens for heavier tails (α = 1.5).
Arcsinh-FM (uniform, with ϕ(x) = arcsinh(x)); TTF and TTFfix (Hickling & Prangle, 2025); gTAF (Jaini et al., 2020). All velocity / coupling networks share the same MLP backbone (4 layers, width 256, SiLU); hyperparameters per dimension were selected by Optuna on a representative configuration (Gumbel, τ = 0.5, α = 2.0) and reused across the grid for that dimension. Metrics. We split metrics by margin type and by joint structure: • Marginal: W1P (mean W1 over Pareto coordinates), W1G (over Gaussian coordinates), Hill estimator on Pareto margins, VaR99 / CVaR99 relative errors, extreme-quantile errors Q99.5 , Q99.9 . • Joint: Absolute Kendall Error (AKE), angular W2 (sliced Wasserstein on the empirical angular measure √ of the top- n extremes), sliced Wasserstein on the full data, and energy distance. We report medians over the 20 replications; “div.” marks runs where the model diverged (W1 > 103 ).
Method
d=10
d=20
d=50
d=100
α = 1.5 Log-FM Arcsinh-FM TTFfix TTF gTAF
0.522 0.561 1.149 0.758 0.534
0.451 0.550 0.963 0.938 0.583
0.534 0.534 0.867 0.758 0.761
0.525 0.592 1.046 0.919 3.320
α = 1.75 Log-FM Arcsinh-FM TTFfix TTF gTAF
0.254 0.262 0.442 0.318 0.270
0.226 0.284 0.417 0.373 0.301
0.269 0.289 0.391 0.354 0.335
0.300 0.316 0.481 0.401 0.934
α = 2.0 Log-FM Arcsinh-FM TTFfix TTF gTAF
0.153 0.166 0.235 0.182 0.171
0.146 0.162 0.222 0.203 0.174
0.195 0.174 0.220 0.207 0.198
0.196 0.212 0.294 0.252 0.219
α = 2.5 Log-FM Arcsinh-FM TTFfix TTF gTAF
0.080 0.082 0.108 0.083 0.088
0.074 0.080 0.107 0.106 0.086
0.090 0.100 0.117 0.113 0.106
0.108 0.116 0.174 0.131 0.225
5.2. Main Results Table 1 reports the headline tail metric, W1P , averaged over the two regular copulas (Gumbel + Gaussian) at the four tail indices and four dimensions. Log-FM is best in 13 out of 16 cells; Arcsinh-FM picks up the remaining 3 at d = 50. Both FM variants are substantially more stable than the baselines, especially for α = 1.5.
5.3. Ablations We isolate three design choices: the Hill gating (adaptive vs uniform), the sampling clamp c, and the ODE solver step count. Hill-gated default vs uniform. Table 4 compares the default Log-FM (with Hill gating, Section 4.2) against the uniform variant and the TTFfix baseline on Gumbel data (d = 20, τ = 0.5, 20 reps). The Hill diagnostic always selects exactly the heavy-tailed coordinates (14 Pareto vs 6 Gaussian out of 20): the threshold αmax = 4 is nonsensitive given α̂ ≈ 1.5–2.5 for Pareto and ≈ 6 for Gaussian. The adaptive variant improves W1G by 30–45% with no W1P cost.
Risk metrics and dependence. Table 2 summarizes the global picture, pooling over both copulas, the three dependence strengths, and all four tail indices. Log-FM is best on every tail-quality and risk metric (WP 1 , CVaR99 , Q99.9 , sliced and energy distances). On the multivariatedependence metrics, TTF/TTFfix retain a small edge at d = 50, 100 while Log-FM is competitive at d = 10, 20. On W1G (Gaussian margins), the Hill-gated default of LogFM closes most of the gap with TTF (further breakdown in Section 5.3).
Clamping (optional). The sampling clamp c is an optional numerical safeguard applied in log-space (corresponding to |x| ≲ ec in data space); it is not used by default. Table 5 reports the relative change in W1P across c values. For every configuration in our main benchmark (α ≥ 1.5) the clamp is inert at c ≥ 10: metrics are bitwise identical to the no-clamp case. The clamp only becomes relevant in extremely pathological tail regimes (α < 1, e.g. in ν = 1/2 Student-t in Appendix E.1), where the population W1 is undefined and exponential error amplification through ϕ−1
Stability. A practitioner cares not only about median performance but about how often the model fails catastrophically. Table 3 reports, by tail index, the fraction of runs with W1P > 1 on Gumbel+Gaussian. Log-FM and ArcsinhFM never exceed 14% even in the most adversarial setting (α = 1.5, d = 10); both reach 0% for α ≥ 1.75. gTAF and TTF suffer divergence rates above 30% in higher dimensions. 8
Tail Annealing for Heavy-Tailed Flow Matching Table 4. Log-FM (uniform ϕ) vs Log-FM (Hill-adaptive) vs TTFfix baseline. Gumbel, d = 20, τ = 0.5, 20 reps, median. Bold = best per α.
Table 2. Median across all Gumbel+Gaussian configurations (480 values per cell). Bold = best per metric per dimension. Method
d=10
d=20
d=50
d=100
W1P (Pareto margins) 0.187 0.216 0.322 0.304 0.304
0.221 0.232 0.308 0.295 0.540
0.233 0.264 0.395 0.349 0.704
CVaR99 (Pareto margins) Log-FM 0.228 0.200 Arcsinh-FM 0.246 0.224 TTFfix 0.880 0.689 TTF 0.494 0.729 gTAF 0.415 0.457
0.257 0.252 0.644 0.687 0.519
0.270 0.373 0.817 0.680 0.705
W1G (Gaussian margins) Log-FM 0.051 0.066 Arcsinh-FM 0.045 0.057 TTFfix 0.033 0.048 TTF 0.035 0.035 gTAF 0.036 0.051
0.080 0.074 0.043 0.030 0.106
0.080 0.085 0.051 0.047 div.
Angular W2 (tail dependence) Log-FM 0.098 0.066 Arcsinh-FM 0.098 0.068 TTFfix 0.148 0.095 TTF 0.121 0.105 gTAF 0.167 0.104
0.071 0.090 0.058 0.046 0.061
0.049 0.049 0.040 0.040 0.098
Log-FM Arcsinh-FM TTFfix TTF gTAF
0.207 0.218 0.322 0.250 0.262
α
Method
W1P
W1G
AKE
1.5
TTFfix Log-FM Log-FM (adaptive)
0.856 0.449 0.468
0.049 0.074 0.043
0.071 0.039 0.044
2.0
TTFfix Log-FM Log-FM (adaptive)
0.225 0.128 0.124
0.034 0.040 0.024
0.066 0.031 0.030
2.5
TTFfix Log-FM Log-FM (adaptive)
0.113 0.068 0.084
0.079 0.110 0.068
0.060 0.030 0.037
Table 5. Clamp ablation: relative W1P change vs. no clamping (c = ∞, default). Gumbel, d = 20, τ = 0.5, 5 reps. The clamp is inert for α ≥ 1.5; it is only relevant in pathological regimes (α < 1).
Table 3. Fraction of runs with W1P > 1 (catastrophic failures), Gumbel + Gaussian, 20 reps. Lower is better. Method
d=10
d=20
d=50
d=100
α = 1.5 Log-FM Arcsinh-FM TTFfix TTF gTAF
14% 11% 58% 32% 8%
2% 6% 50% 45% 24%
4% 2% 36% 36% 42%
2% 2% 54% 46% 59%
α = 2.0 Log-FM Arcsinh-FM TTFfix TTF gTAF
0% 0% 2% 1% 17%
0% 0% 2% 4% 23%
0% 0% 1% 12% 32%
0% 0% 4% 12% 33%
c
α=1.5
α=1.75
α=2.0
α=2.5
5 10, 15, 20 ∞
−2.0% 0.0% 0.0%
−0.7% 0.0% 0.0%
−0.6% 0.0% 0.0%
0.0% 0.0% 0.0%
claim. Dependence is largely preserved by Φ’s diagonal Jacobian; baselines retain only a slight angular-W2 edge at d ≥ 50. Real-data validation on Fama–French (α̂ ≈ 3.7) is in Appendix E.2.
6. Conclusion A coordinate-wise soft-log ϕ(x) = sign(x) log(1 + |x|), alongside a Hill diagnostic, makes standard flow matching work for heavy-tailed data. The mechanism is tail annealing: log-transforming Pareto yields approximately exponential tails, and the induced dynamics implement power transformations X0αt that continuously lighten the tail index. Across 2,880 runs (3 copulas, 4 dimensions, 4 tail indices), Log-FM dominates specialized baselines on W1 , CVaR99 , and extreme-quantile metrics, never diverges severely, and remains competitive on multivariate dependence; on FamaFrench it matches the best baseline. For α < 1 all methods struggle (W1 is undefined); at d ≳ 200 the bottleneck shifts to architecture. Three theoretical directions stand out: a continuous variant via φs2 (Appendix D); combining the transform with discrete normalizing flows, where Q | det JΦ | = j (1 + |xj |)−mj gives exact log-likelihood at O(d) cost; and finite-sample guarantees for log-space score estimation, where bounded scores (Proposition 3.2) place the problem in the regime of existing convergence results.
can produce numerical overflow rather than a meaningful sample. All headline numbers (Tables 1–3) are reported without clamping. ODE solver steps. Sweeping the Euler step count K ∈ {10, 20, 50, 100, 200, 500} on Gumbel (d = 20, τ = 0.5, 10 reps) yields less than 5% variation in W1P ; our default K = 100 is well within the converged regime. The full table is reported in Appendix E.3. 5.4. Discussion The benchmark confirms the theoretical analysis, in that Log-FM has better numerical results whenever marginals are genuinely heavy-tailed (α ≲ 3), with ”few severe divergences in all scenarios” (Table 3) as the strongest practical 9
Tail Annealing for Heavy-Tailed Flow Matching
Impact Statement
Laszkiewicz, M., Lederer, J., and Fischer, A. Marginal tailadaptive normalizing flows. In International Conference on Machine Learning, 2022.
This work targets heavy-tailed generative modeling, with applications in financial risk, insurance, and climate extremes. As with any generative model, outputs should be validated by domain experts before high-stakes use.
Lipman, Y., Chen, R. T., Ben-Hamu, H., Nickel, M., and Le, M. Flow matching for generative modeling. In International Conference on Learning Representations, 2023.
References
Nair, J., Wierman, A., and Zwart, B. The Fundamentals of Heavy Tails: Properties, Emergence, and Estimation. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2022. ISBN 9781316511732.
Albergo, M. S. and Vanden-Eijnden, E. Stochastic interpolants: A unifying framework for flows and diffusions. In International Conference on Machine Learning, 2023. Anderson, B. D. Reverse-time diffusion equation models. Stochastic Processes and their Applications, 12(3):313– 326, 1982.
Pandey, K., Pathak, J., Xu, Y., Mandt, S., Pritchard, M., Vahdat, A., and Mardani, M. Heavy-tailed diffusion models. arXiv preprint arXiv:2410.14171, 2024.
Bengio, Y., Louradour, J., Collobert, R., and Weston, J. Curriculum learning. In International Conference on Machine Learning, pp. 41–48, 2009.
Resnick, S. I. Extreme Values, Regular Variation, and Point Processes. Springer, 1987.
Bingham, N. H., Goldie, C. M., and Teugels, J. L. Regular Variation, volume 27 of Encyclopedia of Mathematics and its Applications. Cambridge University Press, 1987.
Rezende, D. and Mohamed, S. Variational inference with normalizing flows. In International Conference on Machine Learning, pp. 1530–1538, 2015.
Box, G. E. and Cox, D. R. An analysis of transformations. Journal of the Royal Statistical Society: Series B (Methodological), 26(2):211–243, 1964.
Song, J., Meng, C., and Ermon, S. Denoising diffusion implicit models. In International Conference on Learning Representations, 2021a.
Cont, R. Empirical properties of asset returns: Stylized facts and statistical issues. Quantitative Finance, 1(2): 223–236, 2001.
Song, Y. and Ermon, S. Generative modeling by estimating gradients of the data distribution. In Advances in Neural Information Processing Systems, volume 32, 2019.
de Haan, L. and Ferreira, A. Extreme Value Theory: An Introduction. Springer, 2006.
Song, Y., Sohl-Dickstein, J., Kingma, D. P., Kumar, A., Ermon, S., and Poole, B. Score-based generative modeling through stochastic differential equations. In International Conference on Learning Representations, 2021b.
Dinh, L., Sohl-Dickstein, J., and Bengio, S. Density estimation using Real-NVP. In International Conference on Learning Representations, 2017.
Vincent, P. A connection between score matching and denoising autoencoders. Neural Computation, 23:1661– 1674, 2011.
Embrechts, P., Klüppelberg, C., and Mikosch, T. Modelling Extremal Events for Insurance and Finance. Springer Science & Business Media, 2013. Hickling, T. and Prangle, D. Flexible tails for normalizing flows. In Proceedings of the 42nd International Conference on Machine Learning, Proceedings of Machine Learning Research, 2025. arXiv:2406.16971. Ho, J., Jain, A., and Abbeel, P. Denoising diffusion probabilistic models. In Advances in Neural Information Processing Systems, volume 33, pp. 6840–6851, 2020. Hyvärinen, A. Estimation of non-normalized statistical models by score matching. Journal of Machine Learning Research, 6:695–709, 2005. Jaini, P., Kobyzev, I., Yu, Y., and Brubaker, M. Tails of Lipschitz triangular flows. In International Conference on Machine Learning, pp. 4673–4681, 2020. 10
Tail Annealing for Heavy-Tailed Flow Matching
A. Extreme Value Theory: Technical Details This appendix provides the formal definitions and results from extreme value theory used in the main text. See (Resnick, 1987; de Haan & Ferreira, 2006) for a comprehensive treatment. A.1. Heavy-Tailed Distributions Definition A.1 (Heavy-Tailed (Nair et al., 2022)). A random variable X is heavy-tailed if E[eλX ] = ∞ for all λ > 0. Definition A.2 (Regular Variation). A measurable function L : (0, ∞) → (0, ∞) is slowly varying at infinity if L(tx)/L(t) → 1 as t → ∞ for all x > 0. A distribution F is regularly varying with index −α (written F ∈ RV−α ) if its survival function satisfies F̄ (t) = t−α L(t) for some slowly varying L. The Pareto(γ) distribution with F̄ (t) = t−1/γ is regularly varying with index −1/γ (i.e. α = 1/γ); here γ > 0 is the shape parameter (extreme-value index). A.2. The Fisher-Tippett-Gnedenko Theorem Theorem A.3 (Fisher-Tippett-Gnedenko). Let X1 , X2 , . . . be i.i.d. with distribution F . If there exist normalizing sequences an > 0 and bn ∈ R such that maxi≤n Xi − bn P ≤ x → G(x) an for some non-degenerate distribution G, then G must be a generalized extreme value (GEV) distribution: Gξ (x) = exp −(1 + ξx)−1/ξ , 1 + ξx > 0, where ξ ∈ R is the shape parameter (extreme value index). The three cases correspond to: • Fréchet (γ > 0): F has heavy (polynomial) tails; F ∈ RV−1/γ • Gumbel (γ = 0): F has light or sub-exponential tails (e.g. exponential, lognormal) • Weibull (γ < 0): F has finite right endpoint (bounded support). Not to be confused with the Weibull distribution W (k, λ) used in reliability analysis, which has unbounded support and lies in the Gumbel domain; the name clash is historical. Proposition A.4 (Log-transform maps domains). If X is in the Fréchet domain with shape parameter γ > 0, then log X is in the Gumbel domain (γ = 0). This is the qualitative content of Proposition 3.8 and its precise restatement Proposition C.1: regularly varying tails F̄ (t) = t−α L(t) become exponential-type after a logarithm, P(log X > z) = e−αz L(ez ), and distributions with exponential-type tails are in the Gumbel domain (Embrechts et al., 2013, §3.3).
B. Denoising Diffusion Models: Technical Details This appendix presents the framework underlying both denoising diffusion models and flow matching, organized around the stochastic-interpolant formulation of Albergo & Vanden-Eijnden (2023) and Lipman et al. (2023). We then recall how the canonical SDE-based diffusions of Ho et al. (2020) and Song et al. (2021b) arise as the special case in which the interpolation marginals coincide with those of an Ornstein–Uhlenbeck-type forward SDE. B.1. Forward Noising Process The stochastic interpolant of Albergo & Vanden-Eijnden (2023); Lipman et al. (2023) defines the probability path (pt )t∈[0,1] via the marginals pt = Law(Xt ), where Xt = αt X0 + βt X1 ,
X0 ∼ p0 , 11
X1 ∼ N (0, Id ).
(4)
Tail Annealing for Heavy-Tailed Flow Matching
Here X0 and X1 are independent, and (αt )t∈[0,1] and (βt )t∈[0,1] are deterministic schedules such that αt is non-increasing, βt is non-decreasing, with boundary conditions (α0 , β0 ) = (1, 0) and (α1 , β1 ) = (0, 1). This explicit-interpolation viewpoint is not the formalism of Ho et al. (2020); Song et al. (2021b), who instead define the forward process implicitly as the solution to a stochastic differential equation; the two viewpoints are reconciled in Section B.6. From (4), the conditional distribution of Xt given X0 = x0 , denoted qt|0 , is Gaussian: qt|0 (xt | x0 ) = N (xt ; αt x0 , βt2 Id ).
(5)
Common schedules. Two standard choices are the variance-preserving (VP) schedule with αt2 + βt2 = 1 (Ho et al., 2020), which ensures Var(Xt ) = Var(X0 ) when X0 has unit variance, and the linear (flow matching) schedule with (αt , βt ) = (1 − t, t) (Lipman et al., 2023), corresponding to straight-line interpolation between data and noise. Table 6 summarizes common schedule choices with their derivatives. Table 6. Interpolation schedules for flow matching. All satisfy boundary conditions (α0 , β0 ) = (1, 0) and (α1 , β1 ) = (0, 1). The linear schedule corresponds to optimal transport; VP schedules preserve variance when data has unit variance; the quadratic schedule anneals tails more aggressively early in the process. Schedule
αt
βt
α̇t
β̇t
Properties
Linear VP (trig.) VP (poly.) Quadratic
1−t πt cos( √ 2) 1−t (1 − t)2
t πt sin( √2 ) t 1 − (1 − t)2
−1 − π2 sin( πt ) 2 − 2√11−t −2(1 − t)
1 π cos( πt ) 2 2 1 √ 2 t
OT-optimal Smooth endpoints Singular endpoints Fast early annealing
2(1 − t)
B.2. Network Parameterizations: Noise vs. Data Prediction Given a noisy observation xt , two natural conditional-mean targets exist: x̂0 (xt , t) := E[X0 | Xt = xt ],
x̂1 (xt , t) := E[X1 | Xt = xt ].
These correspond to the two standard parameterizations in the diffusion literature (Ho et al., 2020; Song et al., 2021b): (i) Noise prediction (often called the denoiser, or ϵ-prediction in DDPM notation): x̂1 estimates the standard-Gaussian endpoint X1 , which plays the role of the additive noise in (4). (ii) Data prediction (also called the x0 -prediction or x0 -parameterization): x̂0 estimates the clean data sample. The two parameterizations are equivalent in expressive power: given xt = αt X0 + βt X1 from the forward interpolation (4), they satisfy the deterministic identity x̂1 (xt , t) =
xt − αt x̂0 (xt , t) , βt
x̂0 (xt , t) =
xt − βt x̂1 (xt , t) , αt
(6)
so one can be converted into the other at inference time. A third parameterization, v-prediction v(xt , t) := E[α̇t X0 + β̇t X1 | Xt = xt ], is the velocity field used by flow matching and is also a linear combination of x̂0 and x̂1 . B.3. Score-Denoiser Relationship The score of the noised distribution is ∇ log pt (·). For the Gaussian forward kernel (5), the score admits a closed-form expression in terms of the denoiser. Proposition B.1 (Score-Denoiser Identity). Under standard regularity assumptions, h i x̂1 (xt , t) . (7) ∇xt log pt (xt ) = E ∇xt log qt|0 (xt | X0 ) Xt = xt = − βt Proof. The first equality is Fisher’s identity, obtained by exchanging expectation and gradient. For the second, note that ∇xt log qt|0 (xt | x0 ) = (αt x0 − xt )/βt2 . Taking the conditional expectation given Xt = xt : αt X0 − xt αt x̂0 (xt , t) − xt x̂1 (xt , t) E Xt = xt = =− , βt2 βt2 βt where the last equality uses (6). Thus, learning the denoiser x̂1 is equivalent to learning the score ∇ log pt . 12
Tail Annealing for Heavy-Tailed Flow Matching
B.4. Training Objectives The denoiser can be trained by regressing either X0 or X1 from the noised sample Xt = αt X0 + βt X1 . The X0 -prediction loss reads: Z 1 h i 2 LX0 (θ) = EX0 ∼p0 ,X1 ∼N (0,Id ) x̂θ0 (αt X0 + βt X1 , t) − X0 dt, 0
while the X1 -prediction loss is: Z 1 LX1 (θ) =
h
x̂θ1 (αt X0 + βt X1 , t) − X1
EX0 ∼p0 ,X1 ∼N (0,Id )
2
i
dt.
0
Since x̂θ1 = −βt sθ by (7), the X1 -prediction loss is equivalent to denoising score matching (Hyvärinen, 2005; Vincent, 2011). In practice, the integral is approximated via Monte Carlo: sample t ∼ Unif[0, 1], x0 ∼ p0 , x1 ∼ N (0, Id ), form xt = αt x0 + βt x1 , and regress either x0 or x1 . B.5. DDIM Sampling The DDIM framework (Song et al., 2021a) is canonically defined under the variance-preserving constraint αt2 + βt2 = 1 and produces a one-parameter family of reverse transitions sharing the marginals of (5). Given timesteps (tk )K k=0 with tK = 1 and t0 = 0, the transition from tk+1 to tk is: q xtk = αtk x̂θ0 (xtk+1 , tk+1 ) + βt2k − ηt2k x̂θ1 (xtk+1 , tk+1 ) + ηtk z, (8) where z ∼ N (0, Id ) and ηtk ∈ [0, βtk ] so that the radicand is√nonnegative. Translating (Song et al., 2021a) eq. (12) into our √ notation (identifying their σt with ηtk , αt−1 with αtk , and 1 − αt−1 with βtk ) yields exactly (8). The family interpolates between two distinguished members: ηtk = 0 gives the deterministic DDIM update, equivalent to a discretisation of the probability-flow ODE; the choice ηtk = β̃tk , where β̃t2k =
αt2k+1 βt2k 1 − , βt2k+1 αt2k
(9)
is the posterior variance of q(xtk | xtk+1 , x0 ) and recovers ancestral DDPM sampling (Ho et al., 2020; Song et al., 2021a): this matches Song et al. (2021a, eq. (16)) after the same change of variables. The boundary ηtk = βtk is the maximum-noise edge of the family, strictly noisier than DDPM. B.6. Bridging to SDE-Based Diffusion Models The interpolant formulation (4) is not the language used by Ho et al. (2020); Song et al. (2021b), whose forward processes are defined as solutions to discrete- or continuous-time stochastic differential equations. The two formalisms are equivalent at the marginal level for a specific schedule choice, and the equivalence is what Albergo & Vanden-Eijnden (2023) call the “unifying” aspect of stochastic interpolants. We make the bridge explicit in three steps: (i) recover Ho et al.’s DDPM forward chain; (ii) recover Song et al.’s forward score SDE; (iii) recover their reverse-time samplers. (i) DDPM forward chain.
Ho et al. (2020) define a discrete-time Markov chain q DDPM DDPM q(xk | xk−1 ) = N xk ; 1 − βk xk−1 , βk Id ,
with a noise schedule (βkDDPM )K k=1 ⊂ (0, 1). Setting ᾱk =
k = 1, . . . , K,
(10)
DDPM ), marginalising the chain yields j≤k (1 − βj
Q
√ q(xk | x0 ) = N xk ; ᾱk x0 , (1 − ᾱk )Id . √ √ Identifying t = k/K, αt = ᾱk , and βt = 1 − ᾱk recovers (5) exactly. The DDPM noise schedule (βkDDPM ) is therefore one particular choice within the VP family in Table 6, with αt2 + βt2 = 1 by construction. 13
Tail Annealing for Heavy-Tailed Flow Matching
(ii) Forward score SDE.
Song et al. (2021b) unify diffusion models through the continuous-time Itô SDE p dXt = − 12 β(t) Xt dt + β(t) dWt , X0 ∼ p0 ,
(11)
where (Wt )t≥0 is a standard Brownian motion (the “VP-SDE”). The solution at time t has marginals 1 Rt Rt Xt | X0 = x0 ∼ N e− 2 0 β(s) ds x0 , 1 − e− 0 β(s) ds Id , p Rt so identifying αt = exp(− 12 0 β(s) ds) and βt = 1 − αt2 recovers (5). (11) is the continuous-time limit of the DDPM chain (10); both produce the same marginals as (4) under this matched schedule. Albergo & Vanden-Eijnden (2023) show that for any choice of (αt , βt ) in the interpolant, there exists a (possibly time-inhomogeneous) SDE with the same marginals, so the interpolant is strictly more general than the VP family. (iii) Reverse-time samplers. satisfies (Anderson, 1982):
The marginal equivalence carries over to sampling. The time-reversed counterpart of (11) p dX̄t = − 21 β(t)X̄t + β(t)∇x log pt (X̄t ) dt + β(t) dW̄t ,
(12)
with X̄t := X1−t and W̄t a standard Brownian motion. The deterministic probability-flow ODE (Song et al., 2021b) with identical marginals reads dxt = − 12 β(t)xt + 12 β(t)∇x log pt (xt ) dt. (13) DDIM with ηt = 0 is a discretisation of (13); DDPM corresponds to (12). Flow matching learns a velocity field vθ (x, t) ≈ E[α̇t X0 + β̇t X1 | Xt = x] (Lipman et al., 2023); by (7), this velocity is a linear combination of identity and score, so any of the three samplers can be reconstructed from a single trained denoiser regardless of which formalism was used to derive the training loss.
C. Proofs C.1. Precise statement and proof of Proposition 3.8 Proposition 3.8 of the main body states the rate-form conclusion − log P(X̃ > z) = αz + o(z). We give here the corresponding precise two-sided Potter-bound statement and its proof. The result is classical in the regular variation literature; see Bingham et al. (1987, Chapter 1) for a comprehensive treatment. Proposition C.1 (Log-Transform of Regularly Varying; precise). Let X be a nonnegative random variable that is regularly varying with index −α, α > 0. Set X̃ = ϕ(X) with ϕ(x) = sign(x) log(1 + |x|). Then for every ϵ > 0 there exists z0 = z0 (ϵ) such that for all z ≥ z0 , e−(α+ϵ)z ≤ P(X̃ > z) ≤ e−(α−ϵ)z . (14) In particular, − log P(X̃ > z) = αz + o(z) as z → ∞. Proof. By Karamata’s representation, P(X > t) = t−α L(t) for some slowly varying function L. For the soft-log transform ϕ(x) = sign(x) · log(1 + |x|), we have ϕ−1 (z) = ez − 1 for z ≥ 0, and ϕ−1 (z) ∼ ez as z → ∞. Thus P(X̃ > z) = P(X > ez − 1) = (ez − 1)−α L(ez − 1). Since ez − 1 ∼ ez , we have (ez − 1)−α = e−αz (1 + o(1))−α and, by slow variation, L(ez − 1)/L(ez ) → 1. Hence P(X̃ > z) = e−αz L(ez )(1 + o(1)). ′
A standard consequence of slow variation (Resnick, 1987, Proposition 0.8(iv)) is that for every ϵ′ > 0, L(t)/tϵ → 0 and ′ tϵ /L(t) → 0 as t → ∞; equivalently, this follows from Potter’s bounds (Bingham et al., 1987, Theorem 1.5.6) by fixing one argument at a large constant and absorbing the resulting multiplicative factor into the exponent. Either way, for any ϵ′ > 0 there exists t0 = t0 (ϵ′ ) such that for t ≥ t0 , ′
′
t−ϵ ≤ L(t) ≤ tϵ . 14
Tail Annealing for Heavy-Tailed Flow Matching ′
′
Substituting t = ez gives e−ϵ z ≤ L(ez ) ≤ eϵ z for z ≥ log t0 . Combining with P(X̃ > z) = e−αz L(ez )(1 + o(1)) and absorbing the 1 + o(1) factor into a slightly enlarged ϵ (take ϵ′ = ϵ/2 and z0 large enough), we get e−(α+ϵ)z ≤ P(X̃ > z) ≤ e−(α−ϵ)z
for all z ≥ z0 (ϵ),
which is (14). Taking − log and dividing by z yields − log P(X̃ > z)/z → α, i.e. − log P(X̃ > z) = αz + o(z). C.2. Proof of Proposition 3.9 Proof. Let Y = X β . For the survival function of Y : F̄Y (y) = P (X β > y) = P (X > y 1/β ) = F̄X (y 1/β ) Then: F̄X ((ty)1/β ) F̄Y (ty) = F̄Y (t) F̄X (t1/β ) =
F̄X (t1/β · y 1/β ) F̄X (t1/β )
→ (y 1/β )−α = y −α/β as t → ∞, by the regular variation of X with index −α. This extends Lemma 3.4: power transformation preserves the class of regularly varying distributions, with the tail index α scaling inversely with the exponent (equivalently, the shape parameter γ scaling linearly).
D. The Parametrized Soft-Log Family φs2 The main body uses the soft-log ϕ(x) = sign(x) log(1+|x|) uniformly across coordinates, gated by a Hill mask (Section 4.2) when light-tailed margins are present. This appendix develops a continuous one-parameter generalization of ϕ that subsumes both the uniform ϕ and the binary mask as special cases, and clarifies its asymptotic behaviour. The construction is not used in our headline experiments; it is an extension we discuss as a research direction (Section 6). Definition.
For a scale s2 > 0, define φs2 (x) =
1 sign(x) ϕ(s2 x) = log 1 + s2 |x| . s2 s2
The family is parametrized so that φ1 = ϕ (the soft-log of the main body) and φs2 (x) → x as s2 → 0 (the identity). Bulk vs. tail. A Taylor expansion at |x| = 0 gives φs2 (x) = x − 21 s2 x|x| + O(s22 ), so φs2 is the identity to leading order on the bulk |x| ≪ 1/s2 . For |x| ≫ 1/s2 , expanding log(1 + s2 |x|) = log |x| + log s2 + log(1 + 1/(s2 |x|)) gives φs2 (x) =
sign(x) log |x| + log s2 + O (s2 |x|)−1 , s2
which is logarithmic up to an additive constant and an overall 1/s2 scale. Thus s2 acts as the inverse of the cross-over location: the bulk regime extends to |x| ∼ 1/s2 , and the logarithmic-compression regime kicks in beyond. Choosing s2 smaller makes the transform more conservative (closer to identity); choosing it larger compresses more aggressively. Asymptotic tail-mapping. For any fixed s2 > 0, Proposition 3.8 (and its precise restatement Proposition C.1) applies to X̃ = φs2 (X) unchanged up to a constant shift and rescaling: if X ∈ RV−α then − log P(X̃ > z) = α s2 z + o(z)
as z → ∞.
s2 z The proof is identical to that of Proposition C.1 with ϕ−1 (z) = ez − 1 replaced by φ−1 − 1)/s2 . So φs2 still s2 (z) = (e maps regular variation to exponential-type tails for any s2 > 0; only the exponential rate is rescaled by s2 .
15
Tail Annealing for Heavy-Tailed Flow Matching (j)
Coordinate-wise s2 and the adaptive instance. In the multivariate setting we may pick a coordinate-dependent scale (j) (j) s2 . A simple data-driven choice is s2 = c/IQR(Xj ) for a constant c (e.g. c = 1), which puts the cross-over at a robust (j) scale of the marginal. The Hill-gated method of Section 4.2 is the discrete instance of this family with s2 ∈ {0, 1}, (j) (j) chosen by the diagnostic α̂j ≤ αmax ; coordinates flagged heavy-tailed get s2 = 1 (full soft-log) and the others s2 = 0 (j) (identity). An outer scale s1 rescales each X̃j to unit variance, which we apply during preprocessing in both the discrete and continuous variants. Why the binary choice in the main body. Hill-type tail-index estimates are notoriously unstable as point estimates. A (j) continuous s2 derived from α̂j would inject this instability into the preprocessing pipeline; the binary rule uses α̂j only categorically (above or below αmax ), so small perturbations leave the transform unchanged. This matches the stability principle of Section 3.1.
E. Experimental Details This appendix collects the original Hickling Student-t benchmark used at submission time (Section E.1), the real-data Fama–French validation (Section E.2), the NLL validation table against Hickling & Prangle (2025), and the Log-FM hyperparameter listing. The headline benchmark is in Section 5; the experiments below complement it with the original setup the reviewers reference, and are kept for continuity. E.1. Student-t Benchmark iid
We adopt the synthetic benchmark of Hickling & Prangle (2025, Section 4.1): X1 , . . . , Xd−1 ∼ Student-t(ν) and Xd | Xd−1 ∼ N (Xd−1 , 1). We vary d ∈ {10, 20, 50} and ν ∈ {0.5, 1, 1.5, 2, 3, 5, 30} (Table 7 reports ν ∈ {1.5, 2, 3, 5}; the boundary regimes ν ∈ {0.5, 1, 30} are omitted for space). Each configuration uses 5,000 total samples (40/20/40 train/val/test split), averaged over 20 replications. Log-FM achieves the lowest W1 across all configurations. Table 7. Wasserstein-1 distance on the original Hickling benchmark, mean ± standard error over 20 replications. Bold = best. d
ν
TTFfix
TTF
Log-FM
10 10 10 10 20 20 20 50 50 50
1.5 2 3 5 2 3 5 2 3 5
1.78 ± 0.24 0.63 ± 0.10 0.20 ± 0.01 0.12 ± 0.00 0.48 ± 0.03 0.20 ± 0.01 0.13 ± 0.00 0.74 ± 0.04 0.29 ± 0.01 0.17 ± 0.00
17.3 ± 15.2 1.02 ± 0.31 0.90 ± 0.33 0.67 ± 0.20 6.01 ± 4.96 0.95 ± 0.24 1.80 ± 1.34 70.4 ± 64.6 2.73 ± 1.08 14.0 ± 12.7
0.63 ± 0.02 0.25 ± 0.00 0.15 ± 0.00 0.14 ± 0.00 0.23 ± 0.00 0.15 ± 0.00 0.15 ± 0.00 0.28 ± 0.01 0.22 ± 0.01 0.19 ± 0.01
E.2. Real Financial Data We evaluate on the Fama–French 5 industry portfolios (d = 5, daily returns 1963–2023). Hill estimation yields α̂ ≈ 3.7 across dimensions, placing the data in the intermediate-tail regime. Table 8 reports W1 over 10 replications. gTAF achieves a marginal edge; Log-FM is competitive, and substantially better than TTF/TTFfix; mTAF fails catastrophically. On the high-dimensional S&P500 dataset (d = 275) all methods exhibit poor sample quality, dominated by architectural rather than tail-modelling limitations; we omit it. E.3. ODE Solver Step Ablation Table 9 sweeps the Euler step count K at sampling time. Quality plateaus by K = 20 and is essentially flat to K = 500. E.4. Baseline Validation To ensure fair comparison, we validate that our implementations of the baseline methods reproduce the results reported in Hickling & Prangle (2025). Table 10 compares our NLL values (per dimension) with their Table 7 reference values across 16
Tail Annealing for Heavy-Tailed Flow Matching Table 8. W1 on Fama–French 5 (mean ± std, 10 reps). gTAF marginally best; Log-FM competitive; mTAF diverges. Method
W1
gTAF Log-FM TTFfix TTF mTAF
0.127 ± 0.007 0.133 ± 0.013 0.449 ± 0.069 0.487 ± 0.054 5.0 × 104
Table 9. Solver ablation: W1P as a function of Euler step count. Gumbel, d = 20, τ = 0.5, 10 reps. Less than 5% variation between K = 10 and K = 500. α
10
20
50
100
200
500
1.5 2.0 2.5
0.521 0.137 0.084
0.543 0.135 0.081
0.568 0.138 0.083
0.579 0.141 0.084
0.584 0.142 0.085
0.588 0.143 0.086
the original benchmark configurations. Table 10. NLL Validation: Our Implementation vs Hickling Reference (Table 7). Format: ours [ref]. The coupling-layer baselines (TTFfix, TTF) reproduce the reference values within ±0.05 across all configurations; mTAF and gTAF reproduce the reference values closely in moderate-tail regimes (ν ≥ 1) but deviate substantially in pathological regimes (ν = 0.5, especially d = 50), reflecting the training instabilities of these architectures reported in the original paper. d
ν
TTFfix
TTF
mTAF
gTAF
5 5 5 5 10 10 10 10 50 50 50 50
0.5 1.0 2.0 30 0.5 1.0 2.0 30 0.5 1.0 2.0 30
3.32 [3.33] 2.34 [2.35] 1.89 [1.89] 1.47 [1.47] 3.55 [3.54] 2.48 [2.46] 1.94 [1.93] 1.48 [1.47] 3.71 [3.68] 2.58 [2.54] 2.01 [1.98] 1.52 [1.47]
3.32 [3.33] 2.34 [2.34] 1.88 [1.89] 1.47 [1.47] 3.54 [3.55] 2.47 [2.47] 1.94 [1.93] 1.48 [1.47] 3.71 [3.68] 2.58 [2.54] 2.01 [1.98] 1.52 [1.47]
4.05 [4.08] 2.53 [2.49] 1.91 [1.92] 1.46 [1.46] 4.59 [4.48] 2.66 [2.63] 1.96 [1.95] 1.47 [1.46] 6.43 [5.22] 3.14 [2.62] 2.04 [1.98] 1.50 [1.47]
5.70 [6.42] 2.53 [2.49] 1.90 [1.90] 1.47 [1.46] 7.86 [7.13] 2.68 [2.63] 1.95 [1.95] 1.48 [1.47] 21.44 [7.49] 3.66 [2.65] 2.06 [1.99] 1.51 [1.47]
TTFfix and TTF match within ±0.03 across all configurations. mTAF and gTAF show larger deviations in pathological regimes (ν = 0.5, d = 50), consistent with known instabilities reported in the original paper. This validates that our baseline implementations are faithful reproductions. E.5. Hyperparameters Table 11 summarizes all hyperparameters for Log-FM. Baseline hyperparameters follow Hickling & Prangle (2025).
17
Tail Annealing for Heavy-Tailed Flow Matching
Table 11. Log-FM hyperparameters. Network Hidden dimension Layers Activation Time embedding
256 4 SiLU Sinusoidal (dim 256)
Training Optimizer Learning rate Weight decay Batch size Max epochs Early stopping patience Gradient clipping
AdamW 5 × 10−3 10−5 Full batch 5000 100 10.0
Sampling ODE solver Integration steps Output clamp
Euler 100 disabled (c = ∞; optional, inert for α ≥ 1.5)
18