Learning Distributions on Manifolds with Free-Form Flows

Authors: peter  sorrenson
felix  draxler
armand  rousselot
sander  hummerich
ullrich  köthe
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-13 10:14:15

Abstract: We propose Manifold Free-Form Flows (M-FFF), a simple new generative model for data on manifolds. The existing approaches to learning a distribution on arbitrary manifolds are expensive at inference time, since sampling requires solving a differential equation. Our method overcomes this limitation by sampling in a single function evaluation. The key innovation is to optimize a neural network via maximum likelihood on the manifold, possible by adapting the free-form flow framework to Riemannian manifolds. M-FFF is straightforwardly adapted to any manifold with a known projection. It consistently matches or outperforms previous single-step methods specialized to specific manifolds. It is typically two orders of magnitude faster than multi-step methods based on diffusion or flow matching, achieving better likelihoods in several experiments. We provide our code at https://github.com/vislearn/FFF.

Main paper: 
 introduction: 
 Generative models have achieved remarkable success in various domains such as image synthesis [Rombach et al., 2022] , natural language processing [Brown et al., 2020] , scientific applications [Noé et al., 2019] and more. However, the approaches are not directly applicable when dealing with data inherently structured in non-Euclidean spaces, which is common in fields such as the natural sciences, computer vision, and robotics. Examples include earth science data on a sphere, the orientation of real-world objects given as a rotation matrix in SO(3), or data on special geometries modeled by meshes or signed distance functions. Representing such data naively using internal coordinates, such as angles, can lead to topological issues, causing discontinuities or artifacts.Luckily, many generative models can be adapted to handle data on arbitrary manifolds. However, the predominant methods compatible with arbitrary Riemannian manifolds involve solving differential equations-stochastic (SDEs) or ordinary (ODEs)-for sampling and density estimation [Rozen et al., 2021 , Mathieu and Nickel, 2020 , Huang et al., 2022 , De Bortoli et al., 2022 , Chen and Lipman, 2024] . These methods are computationally intensive due to the need for numerous function evaluations during integration, slowing down inference.To address these challenges, we introduce a novel approach for modeling distributions on arbitrary Riemannian manifolds that circumvents the computational burden of previous methods. This is achieved by using a single feed-forward neural network on an embedding space as a generator, with outputs projected to the manifold (Fig. 1 ). We learn this network as a normalizing flow, facilitated by generalizing the free-form flow framework [Draxler et al., 2024 , Sorrenson et al., 2024] to Riemannian manifolds. The core innovation is estimating the gradient of the negative log-likelihood within the tangent space of the manifold.In particular, we make the following contributions:38th Conference on Neural Information Processing Systems (NeurIPS 2024). (Right) We parameterize M-FFF using a neural network in an embedding space, whose outputs are projected to the manifold. This enables simulation-free training and inference, and naturally respects the corresponding geometry, yielding fast sampling and continuous distributions regardless of the manifold.• We extend free-form flows to Riemannian manifolds, yielding manifold free-form flows (M-FFF) in Section 4. • M-FFF can easily be adapted to arbitrary Riemannian manifolds, requiring only a projection function from an embedding space. • It only relies on a single function evaluation during training and sampling, speeding up inference over multi-step methods typically by two orders of magnitude. • M-FFF consistently matches or outperforms previous single-step methods on several benchmarks on spheres, tori, rotation matrices, hyperbolic space and curved surfaces (see Fig. 1 and Section 5). In addition, it is consistently faster than multi-step methods by two orders of magnitude, while also outperforming them in terms of likelihood in several cases.Together, manifold free-form flows offer a novel and efficient approach for learning distributions on manifolds, applicable to any Riemannian manifold with a known embedding and projection. 
related work: 
 Table 1 : Feature comparison of generative models on manifolds. We give a "✓" if any method in a category meets this requirement.Respects topology Single step sampling Arbitrary manifoldsEuclidean ✗ ✓ ✓ Specialized ✓ ✓ ✗ Continuous time ✓ ✗ ✓ M-FFF (ours) ✓ ✓ ✓Existing work on learning distributions on manifolds can be broadly categorized as follows: (i) leveraging Euclidean generative models; (ii) building specialized architectures that respect one particular kind of geometry; and (iii) learning a continuous time process on the manifold. We compare our method to these approaches in Table 1 and give additional detail below.Euclidean generative models. One approach maps the n-dimensional manifold to R n and learns the resulting distribution [Gemici et al., 2016] . Another approach generalizes the reparameterization trick to Lie groups by sampling on the Lie algebra which can be parameterized in Euclidean space [Falorsi et al., 2019] . These approaches come with the downside that a Euclidean representation may not respect the geometry of the manifold sufficiently, e.g. mapping the earth to a plane causes discontinuities at the boundaries. This can be overcome by learning distributions on overlapping charts that together span the full manifold [Kalatzis et al., 2021] . An orthogonal solution is to embed the data and add noise to it in the off-manifold directions, so that the distribution can be learnt directly in an embedding space R m [Brofos et al., 2021] ; this only gives access to an ELBO instead of the exact density. Our method also works in the embedding space so that it respects the geometry of the manifold, but directly optimizes the likelihood on the manifold.Specialized architectures take advantage of the specific geometry of a certain kind of manifold to come up with special coupling blocks for building normalizing flows such as SO(3) [Liu et al., 2023] , SU (d), U (d) [Boyda et al., 2021 , Kanwar et al., 2020] ; hyperbolic space [Bose et al., 2020] ; tori and spheres [Rezende et al., 2020] . Manifold free-form flows are not restricted to one particular manifold, but can be easily applied to any manifold for which an embedding and a projection to the manifold is known. As such, our model is an alternative to all of the above specialized architectures.Continuous time models build a generative model based on parameterizing an ODE or SDE on any Riemannian manifold, meaning that they specify the (stochastic) differential equation in the tangent space [Rozen et al., 2021 , Falorsi, 2021 , Falorsi and Forré, 2020 , Huang et al., 2022 , Mathieu and Nickel, 2020 , De Bortoli et al., 2022 , Chen and Lipman, 2024 , Lou et al., 2020 , Ben-Hamu et al., 2022] . These methods come with the disadvantage that sampling and density evaluation integrates the ODE or SDE, requiring many function evaluations. Our manifold free-form flows do not require repeatedly evaluating the model, a single function call followed by a projection is sufficient.At its core, our method generalizes the recently introduced free-form flow (FFF) framework [Draxler et al., 2024] based on an estimator for the gradient of the change of variables formula [Sorrenson et al., 2024] . We give more details in Section 3.1. 
background: 
 In this section, we provide the background for our method: We present an introduction to free-form flows and Riemannian manifolds.Free-form flows are a class of generative models that generalize normalizing flows to work with arbitrary feed-forward neural network architectures [Draxler et al., 2024] .Euclidean change-of-variables Traditionally, normalizing flows are based on invertible neural networks (INNs, see Kobyzev et al. [2021] for an overview) that learn an invertible transformation z = f θ (x) mapping from data x ∈ R n to latent codes z ∈ R n . This gives an explicitly parameterized probability density p θ (x) via the change-of-variables:EQUATIONwhere f ′ θ (x) ∈ R n×n is the Jacobian matrix of f θ (x) with respect to x, evaluated at x; |f ′ θ (x)| is its absolute determinant. The distribution of latent codes p Z (z) is chosen such that the log-density is easy to evaluate and it is easy to sample from, such as a standard normal. Normalizing flows can be trained by minimizing the negative log-likelihood over the training data distribution:min θ L NLL = min θ E pdata(x) [-log p θ (x)].(2) This is equivalent to minimizing the Kullback-Leibler-divergence between the true data distribution and the parameterized distribution KL(p data ∥p θ ). Sampling from the model is achieved by pushing samples from the latent distribution z ∼ p Z through the inverse of the learned function:x = f -1 θ (z) ∼ p θ .Euclidean gradient estimator Naively computing the volume change log |f ′ θ (x)| in Eq. ( 1) is expensive since it contains the full Jacobian matrix f ′ θ (x) ∈ R n×n , requiring O(n) automatic differentiation steps to compute. Normalizing flow architectures usually avoid this expensive computation by further restricting the architecture to allow fast computation. Luckily, even if such a shortcut is not available, its gradient, which is all we need for training, can still be efficiently estimated as follows:Free-form architectures The central idea of free-form flows is to soften the restriction that the learned model be invertible. Instead, they learn two separate networks, an encoder f θ and a decoder g ϕ coupled by a reconstruction loss, circumventing the need for an invertible neural network f θ :EQUATION)Together, this gives the loss of free-form flows with β, the reconstruction weight as a hyperparameter:EQUATIONThis allows replacing constrained invertible architectures with free-form neural networks. Since f θ is not restricted to efficiently compute the volume change, free-form flows use Eq. ( 5) to compute the gradient of L NLL . To compute Eq. ( 5), free-form flows approximate f -1 θ ′ (z) (which is not tractable) by g ′ ϕ (z) during training:EQUATIONFor density estimation at inference, Draxler et al. [2024] recommend using the explicit decoder Jacobian for the volume change.A manifold is a fundamental concept in mathematics, providing a framework for describing and analyzing spaces that locally resemble Euclidean space, but may have different global structure. For example, a small region on a sphere is similar to the Euclidean plane, but walking in a straight line on the sphere in any direction will return back to the starting point, unlike on a plane.Mathematically, an n-dimensional manifold, denoted as M, is a space where every point has a neighborhood that is topologically equivalent to R n . A Riemannian manifold (M, G) extends the concept of a manifold by adding a Riemannian metric G which introduces a notion of distances and angles. At each point x on the manifold, there is an associated tangent space T x M which is an n-dimensional Euclidean space, characterizing the directions in which you can travel and still stay on the manifold. The metric G acts in this space, defining an inner product between vectors. From this inner product, we can compute the length of paths along the manifold, distances between points as well as volumes (see next section).In this paper, we consider Riemannian manifolds globally embedded into an m-dimensional Euclidean space R m , with n ≤ m. Embedding means that we represent a point on the manifold x ∈ M as a vector in R m confined to an n-dimensional subspace; we write x ∈ M ⊆ R m and denote by π : R m → M a projection from the embedding space to the manifold. A global embedding is a smooth, injective mapping of the entire manifold into R m , its smoothness preserving the topology.In most cases, we work with, but are not limited to, isometrically embedded manifolds, meaning that the metric is inherited from the ambient Euclidean space. Intuitively, this means that the length of a (π(x))) {x ∈ R m : π(x) = x} x → π(x) Rotations SO(d) (d -1)d/2 {Q ∈ R d×d : QQ T = I, det Q = 1} R → arg min Q∈SO(d) ∥Q -R∥ F ; see Eq. (102) Sphere S n n {x ∈ R n+1 : ∥x∥ = 1} x → x/∥x∥ Torus T n = (S 1 ) n n {X ∈ R n×2 : ∥X i ∥ = 1 for i = 1...n} X i → X i /∥X i ∥ for i = 1...n Hyperbolic H n n {x ∈ R n : ∥x∥ < 1} x → x min{1, (1 -ϵ)/∥x∥}path on the manifold is just the length of the path in the embedding space. We note that due to the Nash embedding theorem [Nash, 1956] , every Riemannian manifold has a smooth isometric embedding into Euclidean space of some finite dimension, so in this sense using isometric embeddings is not a limitation. Nevertheless, for some manifolds (especially with negative curvature, e.g. hyperbolic space) there may not be a sensible isometric embedding.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?Answer: [Yes] Justification: The proof for Theorems 2 and 3 are given as Appendices A.1 and A.2 in the appendix. The necessary assumptions are listed with each theorem and we explicitly state in the main text whenever they are relaxed in practice. We intuition on the proofs upon statement in the main text.Guidelines:• The answer NA means that the paper does not include theoretical results.• All the theorems, formulas, and proofs in the paper should be numbered and crossreferenced. • All assumptions should be clearly stated or referenced in the statement of any theorems.• The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition. • Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material. • Theorems and Lemmas that the proof relies upon should be properly referenced. 
manifold free-form flows: 
 The free-form flow (FFF) framework allows training any pair of parameterized encoder f θ (x) and decoder g ϕ (z) as a generative model, see Section 3.1. In this section, we demonstrate how to generalize the steps in Section 3.1 to arbitrary Riemannian manifolds. Note that for simplicity, we choose the same manifold in data and latent spaces, i.e. M X = M Z = M, but the method readily applies to M X ̸ = M Z or G X ̸ = G Z as long as they are topologically compatible, like a sphere and a closed 3D surface without holes. The detailed derivations in the appendix consider this generalization.The volume change on manifolds generalizes the Euclidean variant in Eq. ( 1) by (a) considering the change of volume in the tangent space and (b) accounting for volume change due to changes in the metric: Theorem 2 (Manifold change of variables). Let (M, G) be a n-dimensional Riemannian manifold embedded in R m , i.e., M ⊆ R m . Let p X be a probability distribution on M and let f : M → M be a diffeomorphism. Let p Z be the pushforward of p X under f (i.e., if p X is the probability density of X, then p Z is the probability density of f (X)).Let x ∈ M. Define Q ∈ R m×n as an orthonormal basis for T x M and R ∈ R m×n as an orthonormal basis for T f (x) M.Then, the probability densities p X and p Z are related under the change of variables x → f (x) by the following equation:EQUATION)where Q and R depend on x and f (x), respectively, although this dependency is omitted for brevity.To give an intuition for this result, Fig. 2 shows how the volume change is computed for an isometrically embedded manifold, that is G = I so that |R T GR| = |Q T GQ| = 1. This simplifies Eq. ( 9) to:EQUATION)This is very similar to the familiar change of variables formula in the Euclidean case in Eq. (1), the only difference being that the determinant is evaluated on the n × n projection of f ′ (x) into the tangent spaces. These projections are necessary as the Jacobian of f is singular in the embedding space, since its action is restricted to the local tangent spaces. See the full proof in Appendix A.1.We now generalize the volume change gradient estimator in Theorem 1 to an invertible function on the manifold f θ : M → M. We find that taking the gradient of the manifold change of variables in Eq. ( 9) results in essentially the same computation as in the Euclidean case, but the trace in is now evaluated in the local tangent space: Theorem 3. Under the assumptions of Theorem 2 with f = f θ . Let v ∈ R m be a random variable with zero mean and covariance RR T . Then, the derivative of the change of variables term has the 10) requires to compute the change of a volume element in the tangent spaces under f, which in this example is given by the ratio of lengths of dt and dt ′ . Since f is a map in the embedding space, f ′ (x) defines a linear map between vectors from the embedding space.To correctly compute the change in volume, we use Q and R to change coordinates to the intrinsic tangent spaces, resulting in the linear mapR T f ′ (x)Q : T x M → T f (x) M, which maps dt to dt ′ .following trace expression, where z = f θ (x):EQUATIONThis shows that the adaptation of free-form flows for an invertible function f to isometrically embedded manifolds is remarkably simple (see full proof in Appendix A.2; if the manifold is not isometrically embedded, add the corresponding term in Eq. ( 9)):EQUATIONThe only change is that v must have covariance RR T rather than the identity. We achieve this by sampling standard normal vectors ṽ ∈ R m and then projecting them into the tangent space using the Jacobian of the projection function: v = π ′ (f θ (x))ṽ.(14) Constructing v like this fulfills the conditions of Theorem 3 because E v [v] = 0, and:EQUATION)Just like [Sorrenson et al., 2024 , Draxler et al., 2024] , we further normalize v to reduce the variance of the trace estimator. Equation (13) now allows training invertible architectures on manifolds even if the volume change log |R T f ′ θ (x)Q| is not tractable. Despite using a stochastic estimator for the gradient, we argue in Appendix A.5 that the scaling of the estimator variance with dimension is comparable to the variance due to stochasticity in flow matching and similar methods.Free-form manifold-to-manifold neural networks As discussed in Section 2, invertible architectures have to be specially constructed for each manifold. To overcome this limitation, we now soften the constraint that the learned model be analytically invertible. Instead, we learn a pair of free-form manifold-to-manifold neural networks, an encoder f θ (x) and a decoder g ϕ (z) as arbitrary functions on the manifold:EQUATION)We choose to fulfill Eq. ( 16) using feed-forward neural networks fθ , gϕ : R m → R m working in an embedding space R m of M, but ensure that their outputs lie on the manifold by appending a projection π : R m → M, mapping points from the embedding space R m to the manifold M:EQUATIONFigure 1 illustrates this for an example on a circle M = S 1 .Just like in the Euclidean case, we employ a reconstruction loss to make f θ and g ϕ approximately inverse to one another:L R = E pdata [∥g ϕ (f θ (x)) -x∥ 2 ].(18) We measure the distance in the embedding space; one can modify this to use an on-manifold distance (e.g. great circle distance for the sphere) but we find that ambient Euclidean distance works well in practice, since it is almost identical for small distances and this is the regime we work in. This allows us to substitute f -1 θ ′ (z) ≈ g ′ ϕ (z) in Eq. ( 13):EQUATIONIn Theorem 5 we show that the error of the gradient estimator is bounded by a measure of the mismatch between the encoder and decoder Jacobian matrices. When the encoder and decoder are true inverses, the error reaches zero.Regularization and final loss We find that adding the following two regularizations to the loss improve the stability and performance of our models. Firstly, the reconstruction loss on points sampled uniformly from the data manifold:EQUATIONhelps ensure that we have a globally consistent mapping between the data and latent manifolds in low data regions. Secondly, the squared distance between the output of fθ and its projection to the manifold:EQUATIONdiscourages the output of fθ from entering unprojectable regions, for example the origin when the manifold is S n . The same regularizations can be applied starting from the latent space.The full loss is:L = L M-NLL + β R L R + β U L U + β P L P (22)where the gradient of L M-NLL is computed using Eq. ( 19), and β R , β U and β P are hyperparameters. We give our choices in Appendix B.Question: Does the paper fully disclose all the information needed to reproduce the main experimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?Answer: [Yes] Justification: The gradient estimator in Theorem 3 can be directly implemented using the text below, and we give detailed instructions in Appendix B for each experiment, and released all necessary code. Guidelines:• The answer NA means that the paper does not include experiments.• If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not. • If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.• Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. • While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility.In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. 
experiments: 
 We now demonstrate the practical performance of manifold free-form flows on various manifolds. We choose established experiments to ensure comparability with previous methods, and find:• M-FFF matches or outperforms previous single-step methods. M-FFF uses a simple ResNet architecture, whereas previous methods were specialized to the given manifolds, hindering adoption to novel manifolds. • M-FFF generates samples faster by typically two orders of magnitude than methods sampling in several steps. Despite this great reduction in compute, it achieves a higher generative quality in several cases.In our result tables, we mark as bold (a) the best method overall (both single-and multi-step), and (b) the best single-step method. We provide reconstruction losses of our method and all details necessary to reproduce the experiments in Appendix B. Furthermore, our code is available at https: //github.com/vislearn/FFF. We run each experiment multiple times with different data splits and report the mean and standard deviation of those runs.The group of 3D rotations SO(3) can be represented by rotations matrices with positive determinant, i.e., allQ ∈ R 3×3 with Q T Q = I and det Q = 1.We choose R 3×3 as our embedding space and project to the manifold by solving the constrained Procrustes problem via SVD [Lawrence et al., 2019] (see Appendix B.2).We evaluate M-FFF on synthetic mixture distributions proposed by De Bortoli et al. [2022] with M mixture components for M = 16, 32 and 64. Samples from one of the distributions and samples from our model are depicted in Fig. 3 . Table 3 shows that M-FFF outperforms the normalizing flow developed for SO(3) by Liu et al.[2023], as well as the diffusion-based approaches for the mixtures M = 32 and 64. M = 16 M = 32 M = 64 Fast inference?Moser flow [Rozen et al., 2021] -0.85 ± 0.03 -0.17 ± 0.03 0.49 ± 0.02 ✗: 1000 steps Exp-wrapped SGM [De Bortoli et al., 2022] -0.87 ± 0.04 -0.16 ± 0.03 0.58 ± 0.04 ✗: 500 steps Riemannian SGM [De Bortoli et al., 2022] -0.89 ± 0.03 -0.20 ± 0.03 0.49 ± 0.02 ✗: 100 steps SO(3)-NF [Liu et al., 2023] -0.81 ± 0.01 -0.12 ± 0.004 0.61 ± 0.01 ✓ M-FFF (ours) -0.87 ± 0.02 -0.21 ± 0.02 0.45 ± 0.02 ✓ Earth data on the sphere We evaluate manifold free-form flows on spheres with datasets from the domain of earth sciences. We use four established datasets compiled by Mathieu and Nickel [2020] for density estimation on S 2 : Volcanic eruptions [NGDC/WDS, 2022b], earthquakes [NGDC/WDS, 2022a], floods [Brakenridge, 2017] and wildfires [EOSDIS, 2020] .Figure 1 shows an example for a model trained on flood data. As the reconstruction error sometimes does not drop to a satisfactory level we employ the method described in Appendix B.1 to ensure that the measured likelihoods are accurate. Table 4 shows that M-FFF again outperforms the specialized single-step model; the performance compared to multi-step methods is mixed. We think that multistep models have an advantage on the considered data, as there are large regions of empty space between highly concentrated data points (see density and sample plots in Appendix B.3).Table 4 : M-FFF significantly outperforms the previous single-step density estimator [Peel et al., 2001] on the sphere on real-world earth datasets in terms of negative log-likelihood (lower is better Riemannian CNF [Mathieu and Nickel, 2020 ] -6.05 ± 0.61 0.14 ± 0.23 1.11 ± 0.19 -0.80 ± 0.54 ✗: ∼100 steps Moser flow [Rozen et al., 2021] -4.21 ± 0.17 -0.16 ± 0.06 0.57 ± 0.10 -1.28 ± 0.05 ✗: ∼100 steps Stereographic score-based [De Bortoli et al., 2022] -3.80 ± 0.27 -0.19 ± 0.05 0.59 ± 0.07 -1.28 ± 0.12 ✗: ∼100 steps Riemannian score-based [De Bortoli et al., 2022] -4.92 ± 0.25 -0.19 ± 0.07 0.45 ± 0.17 -1.33 ± 0.06 ✗: ∼100 steps Riemannian diffusion [Huang et al., 2022 To benchmark manifold free-form flows on tori T n , we follow [Huang et al., 2022 ] and evaluate our model on two datasets from structural biology. We consider the torsion (dihedral) angles of the backbone of protein and RNA substructures respectively.We represent a tuple of angles (ϕ 1 , . . . , ϕ n ) ∈ [0, 2π] n by mapping each angle to a position on a circle:X i = (cos ϕ i , sin ϕ i ) ∈ S 1 . Then we stack all X i into a matrix X ∈ R n×2 , compare Table 2.The first dataset is comprised of 500 proteins assembled by [Lovell et al., 2003] and is located on T 2 . The three dimensional arrangement of a protein backbone can be described by the so called Ramachandran angles [Ramachandran et al., 1963] Φ and Ψ, which represent the torsion of the protein backbone around the N -C α and C α -C bonds. The data is split into four distinct subsets General, Glycine, Proline and Pre-Proline, depending on the residue of each substructure.The second dataset is extracted from a subset of RNA structures introduced by Murray et al. [2003] .As the RNA backbone structure can be characterized by seven torsion angles, in this case we are dealing with data on T 7 .We report negative log-likelihoods in Table 5 , finding that M-FFF outperforms a circular spline coupling flow, a normalizing flow particularly developed for data on tori [Rezende et al., 2020] as well as the multi-step methods on one of the datasets. In addition to the quantitative results, we show the log densities of the M-FFF models for the four protein datasets inFig. 5 in Appendix B.4 .We apply M-FFF to the Poincaré ball model, which embeds the 2-dimensional hyperbolic space H 2 of constant negative curvature -1 in the 2-dimensional Euclidean space R 2 , as specified in Table 2 . As this embedding is not isometric, and distances between points grow when moving away from the origin, we must include the last term of Eq. ( 9) when changing variables under a map on this embedded manifold.We show that M-FFF can be applied to non-isometric embeddings using Eq. ( 9) and visualize learned densities in Fig. 1 and in Fig. 6 in Appendix B.5 for several toy datasets defined on the 2-dimensional Poincaré ball model. Further details can be found in Appendix B.5.Manifold with non-trivial curvature Finally, we follow Chen and Lipman [2024] and train M-FFF given by synthetic distributions on the Stanford bunny [Turk and Levoy, 1994] on the data provided with their paper, see Fig. 1 . The natural embedding of this mesh is R 3 , and we train a separate neural network to project from the embedding space to the mesh. This ensures that the projection is continuously differentiable, which we identify to be important for stable gradients.Table 6 shows that M-FFF performs well on this manifold, outperforming Riemannian flow matching in two out of three cases. This experiment underlines the flexibility of our model: We only need a projection function to the manifold in order to train a generative model.Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material? Answer: [Yes] Justification: We provide code to all experiments in Appendix B. Guidelines:• The answer NA means that paper does not include experiments requiring code.• Please see the NeurIPS code and data submission guidelines (https://nips.cc/ public/guides/CodeSubmissionPolicy) for more details. • While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for a new open-source benchmark). • The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https: //nips.cc/public/guides/CodeSubmissionPolicy) for more details. • The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc. • The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only a subset of experiments are reproducible, they should state which ones are omitted from the script and why. • At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable). • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
limitations: 
 Manifold Free-From Flows achieve high generative quality on manifolds despite the approximations made during training: First, the exact inverse of the encoder Jacobian is approximated by the decoder Jacobian, which is implicitly regularized via the reconstruction loss (see Eq. ( 19)). Second, the final gradient computation in Eq. ( 8) is estimated with a single v for each item in the batch, adding noise to the system.At inference time, the negative log-likelihoods we report in all tables are based on the decoder Jacobian. We choose this because even if the decoder ends up not to be invertible after training (that is several latent codes z yield the same generation x = g ϕ (z)), the computed densities are a conservative estimate of the true probability density. The downside is that if the reconstruction loss is high, the likelihoods become inaccurate, see Appendix B.1 for details. We therefore ensure that the final reconstruction losses are vanishing in Table 8 .From a high level perspective, we observe that M-FFF performs less favorable compared to multi-step methods when the density changes sharply or very low density regions are present.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
conclusion: 
 In this paper, we present Manifold Free-Form Flows (M-FFF), a generative model designed for manifold data. To the best of our knowledge, it is the first generative model on manifolds with single-step sampling and density estimation readily applicable to arbitrary Riemannian manifolds. This significantly accelerates inference and allows for deployment on edge devices.M-FFF matches or outperforms single-step architectures specialized to particular manifolds. It also surpasses multi-step methods in several cases, despite reducing the inference compute by typically two orders of magnitude.Adapting M-FFF to new manifolds is straightforward and only requires selecting an embedding space and a projection to the manifold. In contrast, competing multi-step methods are more challenging to adapt as they require implementing a diffusion process or computing distances on the manifold. 
 
 Appendix: In this appendix, we will focus on intuitive definitions of concepts from topology and differential geometry. For a more rigorous treatment of these concepts, see [Jost, 2008] .An n-dimensional manifold M is a space where every point x has a neighborhood which is homeomorphic to an open subset of R n . Intuitively, this means that there is a small region of M containing x which can be bent and stretched in a continuous way to map onto a small region in R n . This is what is meant when we say that the manifold locally resembles R n . If all these maps from M to R n are also differentiable then the manifold itself is differentiable, as long as there is a way to connect up the local neighborhoods in a differentiable and consistent way.The tangent space of the manifold at x, denoted T x M, is an n-dimensional Euclidean space, which is a linearization of the manifold at x: if we zoom in to a very small region around x the manifold looks flat, and this flat Euclidean space is aligned with the tangent space. Because the tangent space is a linearization of the manifold, this is where derivatives on the manifold live, e.g. if f : M X → M Z is a map between two manifolds, then the Jacobianf ′ (x) is a linear map from T x M X to T f (x) M Z .A Riemannian manifold (M, G) is a differentiable manifold which is equipped with a Riemannian metric G : T x M × T x M → R which defines an inner product on the tangent space, which allows us to calculate lengths and angles in this space. The length of a smooth curve γ : [0, 1] → M is given by the integral of the length of its velocity vector γ ′ (t) ∈ T γ(t) M. This ultimately allows us to define a notion of distance on the manifold, as the curve of minimal length connecting two points.In the remainder of the appendix we only consider Riemannian manifolds.Embedded manifolds We define an n-dimensional manifold embedded in R m via a projection functionEQUATIONwhere P ⊆ R m is the projectable set. We require the projection to have the following properties (the first is true of all projections, the others are additional requirements):1. π • π = π 2.π is smooth on P 3. rank(π ′ (π(x))) = n for all x ∈ P Given such a projection, we define a manifold byEQUATIONwith the tangent spaceEQUATIONwhere col denotes the column space. Since the rank of π ′ (x) with x ∈ M is n, the tangent space is n-dimensional and M is an n-dimensional manifold. To avoid clutter we denote the Riemannian metric and its m × m matrix representation with G interchangeably. If M is isometrically embedded then G(x) is just the identity matrix.The Jacobian of the projection is a projection matrix, meaning π ′ (x)π ′ (x) = π ′ (x) for x ∈ M. For any v in the column space of π ′ (x), there is a u such that v = π ′ (x)u and due to the projection property, π ′ (x)v = π ′ (x)u = v. Similarly, for any w in the row space of π ′ (x), wπ ′ (x) = w. If π is an orthogonal projection, π ′ is symmetric by definition and hence the row and column spaces are identical.Integration on embedded manifolds In order to perform integration on the manifold, we cannot work directly in the m-dimensional coordinates of the embedding space, instead we have to introduce some local n-dimensional coordinates. This means that the domain of integration has to be diffeomorphic to an open set in R n . Since this might not be the case for the whole region of integration, we might need to partition it into such regions and perform integration on each individually (each such region, together with its map to R n , is known as a chart and a collection of charts is an atlas).For example, if we want to integrate a function on the sphere, we could split the sphere into two hemispheres and integrate each separately. A hemisphere can be continuously stretched and flattened into a 2-dimensional region, whereas the whole sphere cannot without creating discontinuities.Given an open set U in R n , and a diffeomorphic local embedding function ϕ : U → M, the integral of a scalar function p :EQUATION)The integral on the right is an ordinary integral in R n . The quantity inside the determinant is known as the pullback metric.Theorem 2 (Manifold change of variables). Let (M X , G X ) and (M Z , G Z ) be n-dimensional Riemannian manifolds embedded in R m , i.e., M X , M Z ⊆ R m , and assume they have the same global topological structure. Let p X be a probability distribution on M X and let f : M X → M Z be a diffeomorphism. Let p Z be the pushforward of p X under f (i.e., if p X is the probability density of X, then p Z is the probability density of f (X)).Let x ∈ M X . Define Q ∈ R m×n as an orthonormal basis for T x M X and R ∈ R m×n as an orthonormal basis for T f (x) M Z .Then, the probability densities p X and p Z are related under the change of variables x → f (x) by the following equation:EQUATION)where Q and R depend on x and f (x), respectively, although this dependency is omitted for brevity.Below, we provide two versions of the proof, the second being a less rigorous and more geometric variant of the first.Proof. Let ϕ : R n → M X be defined by ϕ(u) = π X (x + Qu). Let U be an open subset of R n containing the origin which is small enough so that ϕ is bijective. Let ψ : R n → M Z be defined byψ(w) = π Z (f (x) + Rw). Define φ = ψ -1 • f • ϕ and let W = φ(U ). Note that ϕ ′ (u) = π ′ X (x + Qu) • Q and hence ϕ ′ (0) = π ′ X (x)Q = Q (since each column of Q is in T x M X = col(π ′ X (x))). Similarly, ψ ′ (0) = R. Since ψ is a map from n to m dimensions, there is not a unique function from R m to R n which is ψ -1 on the manifold and there are remaining degrees of freedom in the offmanifold behavior which can result in different Jacobians. For our purposes, we define the inverse ψ -1 such that ψ • ψ -1 is an orthogonal projection onto M Z . This means ψ′ (ψ -1 (f (x)))(ψ -1 ) ′ (f (x)) = RR T and hence (ψ -1 ) ′ (f (x)) = R T .Since p Z is the pushforward of of p X under f , the amount of probability mass contained in ϕ(U ) is the same as that contained in f (ϕ(U )) = ψ(W ):EQUATION)and therefore:EQUATION)Changing variables of the RHS with w = φ(u) gives usEQUATION)Since U was arbitrary, we can make it arbitrarily small, demonstrating that the integrands must be equal for u = 0:EQUATIONSince w = ψ -1 (f (ϕ(u))), the Jacobian has the following form when evaluated at the origin (note ϕ(0) = x):EQUATION)Substituting this into the equality, rearranging and taking the logarithm gives the result:EQUATION)Alternative proof Here is a less rigorous and more geometric proof, which may be more intuitive for some readers.Proof. Let x be a point on M X . Consider a small square region U ⊆ M around x (hypercubic region in higher dimensions). If the sides of the square are small enough, the square is approximately tangent to the manifold since the manifold looks very flat if we zoom in. Suppose Q is a basis for the tangent space at x and q 1 , . . . , q n are the columns of Q. Suppose that the sides of the square (or hypercube) are spanned by u i = ϵq i for a small ϵ. The volume spanned by a parallelotope (higher-dimensional analog of a parallelogram) is the square root of the determinant of the Gram matrix of inner products:EQUATIONThe inner product is given by G, namely ⟨u, v⟩ = u T Gv. We can therefore write the volume of U as vol(U ) ≈ ϵ n |Q T GQ|.(36)Now consider how U is transformed under f . It will be mapped to a region f (U ) on M z with approximately straight edges, forming an approximate parallelotope in the tangent space at z = f (x). This region will be spanned by the columns of f ′ (x)ϵQ (since f (x + u i ) ≈ f (x) + f ′ (x)u i ) and hence will have a volume ofEQUATION)where R is a basis for the tangent space at f (x). We can introduce RR T into the expression since it is a projection in the tangent space at f (x) and is essentially the identity within that space. Since the RHS of G Z (z) and the LHS of f ′ (x) both live in this tangent space, we can introduce RR T between them without changing the expression. Then in the last step we use that |AB| = |A||B| for square matrices.The probability density in U and f (U ) should be roughly constant since both regions are very small. Since the probability mass in both regions should be the same we can writeEQUATIONand thereforeEQUATIONwhere the approximation becomes exact by taking the limit of infinitesimally small ϵ. Taking the logarithm, we arrive at the result of Theorem 2.Theorem 3. Under the assumptions of Theorem 2 with f = f θ . Let v ∈ R m be a random variable with zero mean and covariance RR T . Then, the derivative of the change of variables term has the following trace expression, where z = f θ (x):EQUATION)Proof. For brevity, we drop the index θ and denote g = f -1 . First, a reminder that φEQUATION)Note also that since χ(φ(u)) = u, therefore χ ′ (φ(u))φ ′ (u) = I and χ ′ (φ(u)) = φ ′ (u) -1 . Applying Jacobi's formula to φ ′ (u):EQUATIONand substituting in f and g:EQUATIONQ does not depend on θ, but R depends on f (x) and hence θ, so it must be considered in the derivative. However,EQUATIONand since tr(RR T ) = tr(I) is a constant, tr(R∇ θ R T ) = 0. Expanding Eq. ( 47):∇ θ log |R T f ′ (x)Q| = tr(∇ θ (R T )f ′ (x)QQ T g ′ (f (x))R) + tr(R T ∇ θ (f ′ (x))QQ T g ′ (f (x))R).(49) Since Q is an orthonormal basis for T x M X , QQ T is a projection matrix onto T x M X . This is because (QQ T ) 2 = QQ T QQ T = QQ T , using Q T Q = I. As a result, QQ T π ′ (x) = π ′ (x). Since g can also be written inside a projection: g(z) = π Z (g(z)), therefore g ′ (z) = π ′ Z (g(z))g ′ (z), so QQ T g ′ (z) = g ′ (z). Note also that f ′ (x)g ′ (f (x)) = I since f • g = id. This simplifies the equation:EQUATIONand finallyEQUATIONIn the above proof we used the fact that QQ T g ′ (z) = g ′ (z), where we dropped the index θ and use g := f -1 for brevity. Can we use RR T f ′ (x) = f ′ (x) to simplify the equation further? No, we cannot, since the expression involving f ′ is actually its derivative with respect to parameters, which may not have the same matrix structure as f ′ . Is it instead possible to use g ′ (z)RR T = g ′ (z) for simplification? If g is implemented as π Z (g(z)), this is not necessarily true, as g ′ (z) might not be a map from the tangent space at z to the tangent space at g(z). For example, if we add a small deviation v to z, where v is orthogonal to the tangent space at z, then g(z + v) might not equal g(z). However, this would mean that derivatives in the off-manifold direction can be non-zero, meaning that g ′ (z)v ̸ = g ′ (z)RR T v = 0 (since RR T will project v to 0). We can change this by prepending g by a projection:EQUATION)If π Z is an orthogonal projection, meaning that π ′ Z is symmetric, the column space and row space of π Z will both be the same as those of RR T , meaning π ′ Z (z)RR T = π ′ Z and hence g ′ (z)RR T = g ′ (z). This leads to the following corollary: Corollary 4. Suppose the assumptions of Theorem 2 hold with f = f θ and the following implementation:f -1 θ = π X • f -1 θ • π Z (53)where π Z is an orthogonal projection. Then the derivative of the change of variables term has the following trace expression, where z = f θ (x):EQUATION)Proof. Again, we drop the index θ and let g = f -1 for brevity. Take the result of Theorem 3 and use the cyclic property of the trace and the properties of g ′ discussed above:EQUATION)We use Hutchinson-style trace estimators to approximate the traces given above. This uses the property that, for a matrix A ∈ R n×n and a distribution p(v) in R n with unit second moment (meaningEQUATIONmeaning that v T Av ≈ tr(A) is an unbiased estimate of the trace of A.We have two variants of the trace estimate derived above, one evaluated in R n , the other in R m . The first can be estimated using the following equality (again dropping the index θ and using g = f -1 ):EQUATION)where p(u) has unit second moment in R n and p(v) is the distribution of Ru, which lies in the tangent space at x and has unit second moment in that space by which we mean E[vv T ] = RR T . An example of such a distribution is the standard normal projected to the tangent space, i.e. v = RR T ṽ where ṽ is standard normal.In the second case, we can just sample from a distribution with unit second moment in the embedding space R m :EQUATIONThe error bound on the gradient from Draxler et al. [2024, Theorem 4.2 ] can be readily extended to Riemannian manifolds: Theorem 5. Carry over the assumptions of Theorem 2 with f = f θ and let g ϕ be a manifold-tomanifold function. LetJ f θ = f ′ θ (x), J g ϕ = g ′ ϕ (z), and J f -1 θ = f -1′ θ (z). Then: tr(R T (∇ θ J f θ )J g ϕ R) -∇ θ log |R T J f θ Q| ≤ ∥R T (∇ θ J f θ )J f -1 θ R∥ F ∥R T J f θ J g ϕ R -I n ∥ F . (64)Proof. The proof closely follows [Draxler et al., 2024] and utilizes the Cauchy-Schwarz inequality for the Frobenius inner product, which states that for matrices A and B, we have | tr(A T B)| ≤ ∥A∥ F ∥B∥ F . Applying this to our case:EQUATION)We can re-express this term by introducing the identity matrix in terms of the Jacobians of f θ and its inverse:EQUATION)By applying the Cauchy-Schwarz inequality, we obtain the bound:EQUATION)To further clarify, we recall the function φ introduced in the proof of Theorem 2. This yields:EQUATION)Thus, we can represent J g ϕ as:EQUATIONwhere we used J g ϕ = QQ T J g ϕ and J f -1 θ = QQ T J f -1 θ using similar reasoning to the proof of Theorem 3.When using a Hutchinson trace estimator with standard normal v ∈ R n , we can reduce the variance of the estimate by scaling v to have length √ n (see [Girard, 1989] ). The scaled variable will still have zero mean and unit covariance so the estimate remains unbiased, but the variance is reduced, with the effect especially pronounced in low dimensions.While we can take advantage of this effect in both our options for trace estimator, the effect is more pronounced in lower dimensions, so we reduce the variance more by estimating the trace in an n-dimensional space rather than an m-dimensional space. Hence the first version of the trace estimator, where v is sampled from a distribution in T x M X is preferable in this regard.Let's provide some intuition with an example. Suppose n = 1, m = 2 and R = (1, 0) T . We want to estimate the trace of A = diag(1, 0). Using the first estimator, we first sample v = RR T ṽ with ṽ standard normal which results in v = (u, 0) T where u ∈ R is standard normal. Then we scale v so it has length √ n = 1. This results in v = (r, 0) T where r is a Rademacher variable (taking the value -1 and 1 with equal probability). The trace estimate is therefore r 2 = 1, meaning we always get the correct answer, so the variance is zero. The second estimator samples v directly from a 2d standard normal, then scales it to have length √ m = √ 2. Hence v is sampled uniformly from the circle with radius √ 2. We can write v = √ 2(cos θ, sin θ) T with θ sampled uniformly in [0, 2π]. The estimate v T Av = 2 cos 2 θ. This is a random variable whose mean is indeed 1 as required but has a nonzero variance, showing that the variance is higher when estimating in the m-dimensional space.For this reason, we choose the first estimator, sampling v in the tangent space at x. This also simplifies the definition of g, meaning that we don't have to prepend it with a projection.Since generating in high-dimensional spaces can raise concerns about an estimator's scaling behavior, we argue theoretically that the free-form flow estimator scales comparably to flow matching as dimensions increase. Experimental results further support this, showing that non-Riemannian freeform flows exhibit strong scaling performance in spaces up to 261 dimensions [Draxler et al., 2024] .One metric to assess scaling behavior is to consider the variance of the gradient estimator. With standard normal noise, the variance of the trace estimator v T Av is 2∥A∥ 2 F [Hutchinson, 1989] . If we apply this to a simple linear free-form flow model f (x) = Ax and g(z) = Bz, we find that summing up the variances of the gradient estimates with respect to each element of A leads to a total variance of 2n tr(BB T ). Note that in a converged model, BB T is equal to the covariance of the data. A similar calculation for a flow matching loss of the form 1 2 ∥Ax -y∥ 2 leads to a total variance ∥x∥ 2 tr(Σ)with Σ the covariance of p(y|x). These results are proven below. We can see that both expressions scale as n 2 , assuming that ∥x∥ 2 and the trace terms scale as n. These assumptions are fulfilled if, for example, the data is Gaussian, with covariance that doesn't depend on n. Since the number of parameters (elements of A) scales as n 2 , the variance per parameter is constant. We thus expect similar scaling behavior to flow matching, with no problems due to variance in high dimensions. Lemma 6. Let A ∈ R n×n be a matrix with entries A ij , and let B ∈ R n×n be any matrix. DefineC ij = ∂A ∂A ijB for each i, j. Let v ∈ R n be a random vector with entries independently drawn from the standard normal distribution. Then, the total variance of the Hutchinson estimators v T C ij v for tr(C ij ) over all i, j is given byEQUATIONwhere ∥ • ∥ F denotes the Frobenius norm.Proof. We know from Hutchinson [1989] thatVar v T C ij v = 2∥C ij ∥ 2 F (75)Note the form of ∂A ∂Aij , using the Kronecker delta:EQUATIONThis impliesEQUATIONThus, C ij has non-zero entries only in the i-th row, and that row is equal to the j-th row of B:EQUATION)Next, calculate the Frobenius norm of C ij :EQUATIONwhere B j: denotes the j-th row of B.Finally, sum over all i and j to find the total variance:EQUATIONSinceEQUATIONthe total variance simplifies to Total Variance = 2n∥B∥ 2 F .(82) Noting that ∥B∥ 2 F = tr(BB T ) completes the proof. Lemma 7. Let A ∈ R n×n be a fixed matrix, x ∈ R n a fixed vector, and let y ∈ R n be a random vector with conditional distribution p(y|x) having covariance matrix Σ. Define the loss functionL(A) = 1 2 ∥Ax -y∥ 2 .Then, the total variance of the gradient estimators ∂L ∂A ij over all i and j under p(y|x) is given by Total Variance = ∥x∥ 2 tr(Σ).Proof. The loss function is given byEQUATIONThe derivative of L(A) with respect to A ij isEQUATION)Since Ax is a vector whose k-th component is (AxEQUATIONwhere δ ki is the Kronecker delta.Therefore, the derivative becomesEQUATIONWe are interested in the variance of the estimator ∂L ∂A ij under p(y|x). Since A and x are fixed, the only randomness comes from y. Assuming that y is a random vector with mean µ = E[y|x] and covariance matrix Σ, we haveEQUATIONSince (Ax -y) i = (Ax) i -y i , and (Ax) i is deterministic, it follows thatEQUATIONTherefore,EQUATIONThe total variance over all i and j isEQUATIONTherefore, we have proven that the total variance is ∥x∥ 2 tr(Σ).In accordance with the details provided in Section 4, our approach incorporates multiple regularization loss components in addition to the negative log-likelihood objective. This results in the final loss expression:EQUATIONFor each of the terms, there is a variant in xand in z-space, as indicated by the superscript. In detail:The first loss L R represents the reconstruction loss:EQUATIONEQUATIONHere d(x, y) = ∥x -y∥ 2 is the standard reconstruction loss in the embedding space. This could be replaced with a distance on the manifold. However, this would be more expensive to compute and since we initialize networks close to identity, the distance in the embedding space is almost equal to shortest path length on the manifold.In order to have accurate reconstructions outside of training data, we also add reconstruction losses for data uniformly sampled on the manifold, both for x and z:EQUATIONEQUATIONFinally, we make sure that the function learned by the neural networks is easy to project by regularizing the distance between the raw outputs by the neural networks in the embedding space and the subsequent projection to the manifold (compare Eq. ( 17)):EQUATIONEQUATIONIf these superscripts are not specified explicitly in the following summary of experimental details, we meanβ x = β z .In all cases, we selected hyperparameter using the performance on the validation data.Sampling from a trained model can be easily achieved by sampling from the latent distribution and performing a single pass through the decoder g. However, in order to evaluate the likelihood of the test set under our model, as in [Draxler et al., 2024] , we need to calculate the change of variable formula w.r.EQUATIONEQUATION)Here we used the approximation f ≈ g -1 . While is expensive to compute for training, it is reasonably fast to compute during inference time. To show the validity of this approximation, we report the reconstruction losses computed on the test dataset for all experiments in table 8.Insufficient reconstruction losses. If the assumption of f ≈ g -1 is not sufficiently fulfilled the measured likelihoods might be inaccurate. We try to identify such cases by sampling points from M around the proposed latent point with small noise strength σEQUATIONAs inverse of the decoder g we use the sample z which results in the lowest reconstruction lossEQUATIONWe also do this with samples drawn around x = π(x + N (0, σ)) and uniformly from the manifold z = U(M). We note that in most cases except for the earth datasets the likelihoods computed in this way agree with f ≈ g -1 . In case of the earth datasets, we note that the newly computed likelihoods now agree with observed model quality. Specifically, whenever reconstruction loss is low we also see agreement between the likelihoods computed via f and our approximation of g -1 via sampling. Otherwise, the disagreement between f and g -1 and the resulting drop in model quality are correctly diagnosed. Therefore, for the all S 2 experiments we report the likelihoods computed by our approximation.Sampling based evaluation. In order to evaluate our models independently of their reconstruction capabilities, we propose to also report a sample-based metric. Due to the curse of dimensionality, sample-based metrics are only tractable in low dimensions. Therefore, we report the Wasserstein-2 distance between model samples and the test dataset, including the standard deviation over multiple training runs, for all 2-dimensional manifolds in Table 7 . As competing method have neither reported a sampling based evaluation metric nor published their models, we propose to use our baseline results for future benchmarking. This justifies evaluating M-FFF via negative log-likelihoods.SO(3) (K = 16) (2.7 ± 0.3) × 10 -5 SO(3) (K = 32) (1.6 ± 1.8) × 10 -4 SO(3) (K = 64) (7.2 ± 1.6) × 10 -5(5.1 ± 0.8) × 10 -6 Earthquakes (1.9 ± 0.8) × 10 -7 Flood (1.9 ± 0.6) × 10 -5 Fire (1.7 ± 0.4) × 10 -6 General (7.2 ± 0.9) × 10 -6 Glycine (1.6 ± 0.2) × 10 -5 Proline (1.3 ± 0.3) × 10 -5 Pre-Pro(1.5 ± 0.8) × 10 -4 RNA (8.9 ± 0.6) × 10 -4Bunny (k = 10) (2.6 ± 0.5) × 10 -5 Bunny (k = 50) (2.2 ± 0.5) × 10 -5 Bunny (k = 100) (2.3 ± 0.5) × 10 To apply manifold free-form flows, we project an output matrix R ∈ R 3×3 from the encoder/decoder to the subspace of special orthogonal matrices by finding the matrix Q ∈ SO(3) minimizing the Frobenius norm ∥Q -R∥ F . This is known as the constrained Procrustes problem and the solution Q can be determined via the SVD R = U ΣV T [Lawrence et al., 2019] :EQUATIONwhere the diagonal entries of Σ were sorted from largest to smallest and S = Diag(1, . . . , 1, det(U V T )).The special orthogonal group data set is synthetically generated. We refer to [De Bortoli et al., 2022] for a description of the data generation process. They use an infinite stream of samples. To emulate this, we generate a data set of 10 7 samples from their code, of which we reserve 1,000 for validation during training and 5,000 for testing. We vectorize the 3 × 3 matrices before passing them into the fully-connected networks. All training details are given in Table 9 , one training run takes approximately 7 hours on a NVIDIA A40.The data set is synthetically generated; of the N = 100, 000 data points, we use 1% of for validation and hyperparameter selection and 5% for test NLL. Each run uses a different random initialization of papers.We follow previous works and use a dataset split of 80% for training, 10% for validation and 10% for testing. For the earth datasets we use a mixture of 5 learnable Von-Mises-Fisher distributions for the target latent distribution. We base our implementation on the hyperspherical_vae library [Davidson et al., 2018] . In order to stabilize training we apply a small amount of Gaussian noise to 12 : Details on the torus datasets. Each dataset is randomly split into a train dataset (80%), validation dataset (10%) and test dataset (10%). During training, we add Gaussian noise with mean zero and standard deviation given by 'noise strength' to the data, to counteract overfitting. every batch (see table 10 ) and project the resulting data point back onto the sphere. Other training hyperparameters can be found in table 11. Each model trained around 20h on a compute cluster using a single NVIDIA A40.The torus datasets are randomly split into a train dataset (80%), validation dataset (10%) and test dataset (10%). To counteract overfitting, we augment the RNA dataset with random Gaussian noise. The noise strength and total number of instances is reported in 13 : Details on the model architecture, loss weights and optimizer parameters for the torus datasets. We use the same configuration for all protein datasets on T 2 . Figure 5 : Log density of M-FFF models in the (Φ, Ψ)-plane of protein backbone dihedral angles (known as Ramachandran plot [Ramachandran et al., 1963] ). The learned density matches the true density indicated by the test dataset (black dots) very well. Note also that the learned distribution obeys the periodic boundary conditions. distribution. We train for 120k steps with a batch size of 512 which takes 2.5 to 3 hours on a NVIDIA GeForce RTX 2070 graphics card. Further hyperparameters used in training can be found in Table 13 .A straightforward way to define distributions on hyperbolic space (but also other Riemannian manifolds) is, to define a probability density p tangent in the tangent space at the origin and use the exponential map exp 0 to pushforward this distribution onto the manifold using Eq. ( 9):EQUATIONwhere G manifold denotes the metric tensor of the embedded manifold and G tangent the metric tensor of the tangent space. This is also known as a 'wrapped' distribution.We use the Poincaré ball model, which embeds the n-dimensional hyperbolic space H n in the ndimensional Euclidean space R n as defined in Table 2 . The exponential map at the origin of this embedding and its Jacobian determinant are simply given by:EQUATION)The metric tensor at some point p ∈ H is defined by: G ij H (p) = λ 2 p δ ij with λ p = 2/(1 -∥x∥ 2 ). The metric tensor of the tangent space is the usual Euclidean metric tensor: G ij R = δ ij . With this at hand, we can define latent and toy distributions as wrapped distributions at the origin, as depicted in Fig. 6 .For training, we sample 100k data points from each distribution. Hyperparameters for each model can be found in Table 14 . Training takes approximately 2.5, 16, 8 and 16 hours on a NVIDIA A40 graphics card for the one Gaussian, five Gaussians, swish and checkerboard dataset respectively. The resulting model densities are shown in Fig. 6 .We base our experiment on the manifold and data provided by [Chen and Lipman, 2024] using 80% for training, 10% for validation and hyperparameter tuning. We report test NLL on the remaining 10% of the data. Each run starts from different parameter initialization. They give the manifold as a triangular mesh, consisting of vertices v i ∈ R 3 , i = 1, . . . N v and triangular faces f j ∈ {1, . . . , N v } 3 .Since the projection to the nearest point on the mesh has zero gradient in parts of R 3 , we instead project to the manifold using a separately trained auto encoder with a spherical latent space. This autoencoder consists of an encoder e : R 3 → R 3 consisting of five hidden layers with 256 neurons each, SiLU activations and an overall skip connection. The latent codes are computed by projecting the encoder outputs e(x) to a sphere as z(x) = e(x)/∥e(x)∥, so that the latent space has the same the topology as the input mesh. Then, a decoder d(z) with the same structure as the encoder is trained to reconstruct the original points by minimizing ∥x -d(e(x)/∥e(x)∥∥ 2 . We train it for 2 18 = 262, 144 steps, with each batch consisting of all N v = 2, 502 vertices and an additional N = 2, 502 uniformly • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 


Decision: The provides a simple extension of free-form flows to manifolds. While the work is incremental (the driving theory was available in the recent papers by Draxler & Sorrenson) it nonetheless makes the ideas explicit and the empirical demonstrations are valuable. The authors are strongly encouraged to update the paper according to the key discussions, in particular, the low reconstruction loss was generally deemed an important concern. The main empirical observation is the improvement of single-step methods; it would perhaps be reasonable to tone down the claims regarding the comparison to multi-step methods as the submitted text could be misunderstood to imply that the proposed method is competitive with such.

Reviews: review: paper_summary: paper_summary: This paper introduces Manifold Free-Form Flows (M-FFF), an extension of free-form flows to Riemannian manifolds. M-FFF is easy to adapt to an arbitrary Riemannian manifold, since the method requires only a projection function from an embedding space. Moreover, the method itself relies only on a single function evaluation during training and sampling, which speeds up inference over multi-step methods by several orders of magnitude. As demonstrated by the experimental section, M-FFF consistently matches or outperforms single-step methods on several benchmarks on spheres, tori, rotation matrices, and hyperbolic space, while being reasonably competitive with multi-step methods.
strength_weakness: strength_weakness: 1. The method, in essence, trades the hard manifold constraints and theoretical cleanliness of previous Riemannian Normalizing Flows (Mathieu et al., 2020) for a projection-based approximation. Specifically, the benefit of a loose and adaptable "encoder-decoder" framework comes at the cost of having to optimize three additional losses: the reconstruction loss $\mathcal{L}_R$, the uniformity loss $\mathcal{L}_U$, and the projection loss $\mathcal{L}_P$. That is to say, unless these losses are optimized well, the geometry may be respected very approximately, which is a drawback relative to existing multi-step methods. The authors themselves admit that they have had issues with the reconstruction loss being too high (see Section 5). Additionally, such a method may be rather suboptimal for applications where there are hard physical constraints, e.g. see learning symmetric densities over $SU(n)$ [a] in the context of lattice QFT for downstream use in gauge integrals. M-FFF is nonetheless an interesting method, but these limitations should be well-noted.

2. There are clearly contexts in which this method is suboptimal relative to existing multi-step methods, e.g. see Figure 6. Due to the many approximations made in the method, it can be hard to match some of the more rigid density geometry exhibited by the checkerboard pattern, that strict on-manifold multi-step methods have no trouble learning (see Lou et al. (2020)). This should also be well noted; the speed flexibility of M-FFF come at a cost.

### References

[a] Equivariant manifold flows. https://arxiv.org/abs/2107.08596
questions: questions: My questions to the authors are listed below:

1. In Section 5, the authors mention that the "reconstruction error sometimes does not drop to a satisfactory level." My question is, when does this usually happen, and can we more precisely characterize the situations in which one of the non-NLL losses fails to give a satisfactory result? I believe this question is key to better determining when precisely M-FFF is a reasonable and/or good alternative to strict multi-step methods such as Lou et al. (2020) or Mathieu et al. (2020).

### Additional Comments and Minor Corrections

The writing in this paper is for the most part quite clear and the presentation of material is excellent, with several expository figures given to elucidate important parts of the method. Nonetheless, I would like to give a non-exhaustive list of minor corrections below:

L62: "for building normalizing flows such as SO(3)" -> "for building normalizing flows for manifolds such as SO(3)"

L203: I believe the authors meant to say the substitution happens for equation (13), and the results is shown in equation (19). That is, "This allows us... in Eq. (19)" -> "This allows us.... in Eq. (13)"

L228: "by rotations matrices with positive determinant" -> "by rotation matrices with positive determinant"
limitations: limitations: Although some limitations of the method are mentioned throughout the paper, this work would benefit from an explicit section made to address limitations of the method (in particular incorporating some of what I have mentioned in the "Weaknesses" section above).
paper_summary: paper_summary: The extends a recent class of generative models, called free-form flows, to Riemannian manifolds.  These models are trained like normalizing flows, but relax the bijective constraint and instead employ separate encoder and decoder networks that are trained to be inverses of each other.  The proposed method differs from existing flows for manifolds in that it constructs a map from the manifold M, back to M, using a full rank transformation plus a projection.  This differs from existing manifold flows which typically use injective transformations.  The authors derive a change of variables on manifolds for this setting and derive an unbiased estimator of its gradient and propose a loss function that affords stable training of these models.  The method consistently outperforms existing generative models on manifolds at learning the various target distributions from samples at a significantly lower expense than previous methods.
strength_weakness: strength_weakness: - I'm not fully convinced by the likelihood based evaluation.  The invertibility of normalizing flows ensures that the flow can produce samples anywhere on $\mathcal{M}$ and similarly, can assign a likelihood to any point on $\mathcal{M}$.  Since free form flows are not invertible, it is possible that the flow does not map to all of $\mathcal{M}$.  This would make the comparison between the likelihood values of normalizing flows and free form flows invalid because the likelihoods are valid on different domains, which could potentially invalidate the reported results.  Instead, the experiments could report a sample based metric to measure the method's performance.
- There is no error bound on the likelihood gradient like there is in Draxler et al. 2024.  Although it seems like the method performs well in practice, I would like to know exactly how much is lost by using a decoder network instead of a true inverse.  In injective flows, I don't think that the free form flow gradient works because $\nabla \log |f'(x)| = \text{Tr}{(\nabla f'(x){f^+}'(z))}$ where ${f^+}'$ is the pseudo inverse of $f'$, so using any other possible inverse would give a wrong gradient.  I understand that the proposed method is not an injective flow and might not suffer from this, but I would like to see a gradient error bound to know for sure.
- The $m$ and $n$ dimensions in the proof of theorem 2 seem to be flipped in the appendix.
questions: questions: - It would be nice to clarify in the text that all of the quantities considered are written in terms of the Euclidean ambient space coordinates.
limitations: limitations: The authors adequately addressed the limitations.
paper_summary: paper_summary: This paper introduces M-FFF, a single-step generative model that operates on Riemann manifolds. The authors extend the existing FFF, which functions in Euclidean space, establishing the theoretical foundation for learning distributions on various Riemann manifolds and enabling faster generation than existing multi-step and single-step methods. The effectiveness of the proposed approach is evaluated on several datasets.
strength_weakness: strength_weakness: 1. Extending existing generative models operating in Euclidean space to Riemannian manifold is a well-trodden path in recent research, so this extension alone lacks novelty. However, the theoretical contributions specific to the extension are solid.
2. The authors evaluated the method's effectiveness only under limited conditions where the dimensionality n is ~10 or less. Given that the theory and method are highly general without specific applications, there should be a discussion on the potential limitations of the algorithm.
questions: questions: 1. Is it realistic to generate data in high-dimensional spaces (e.g., 100 ~ 1000 dimensions), or is it inherently challenging? Understanding the short rebuttal period, I don't expect new experiments, but could you provide insights on the expected impacts of increased dimensionality?
2. Do you plan to release the source code, including reproduction experiments? Considering the scarcity of source code for generative models on Riemann manifolds, releasing the code would significantly contribute to the field.
3. Regarding Specialized architectures in related work, refer to, for example, the paper for quotient manifolds by Yataka et al., Grassmann manifold flows for stable shape generation, NeurIPS2023. While authors of M-FFF claim their theoretical results operate on arbitrary Riemann manifolds, can the same be said for quotient manifolds?
4. M-FFF is described as requiring "requiring only a projection function from an embedding space." However, isn't the calculation of R and Q dependent on the manifold necessary for computing the negative log-likelihood? If so, how are these obtained? At least, it seems there is a discrepancy between the claims in the Introduction and the necessity of Eq. (9) in Section 5 (Toy distributions on hyperbolic space).
5. The authors describe that evaluating the loss in Euclidean space (Eq. (22)) is practically acceptable because "we find that ambient Euclidean distance works well in." However, considering that Riemannian flow matching is rigorously formulated based on the Riemannian metric, this approach seems inconsistent. Theoretically, wouldn’t a formulation using the Riemannian manifold's metric be necessary?
6. Particularly in spaces with negative curvature, such as hyperbolic space, I am concerned about the appropriateness of using Euclidean distance. If some measures are needed, please explain why the formulation in Eq. (22) is appropriate as a general evaluation metric for M-FFF.
limitations: limitations: 1. There is no mention of limitations related to dimensionality.
paper_summary: paper_summary: In this paper the authors propose to adapt the formalism in Draxtler et al. (2024) and Sorrensen et al. (2024) to Riemannian geometry. The idea they build upon is called Free-form flows and consists of: (i) relaxing the bijectivity requirement of Normalizing Flows by approximating the inverse transformation with a learnable function and (ii) approximating the gradient of the Jacobian determinant with the Hutchinson's trace estimator (and in particular to use the approximate inverse at this point). The contribution consists in generalizing such framework to accommodate for Riemannian manifolds (and specifically to accommodate for the further multiplication by the tangent space orthonormal basis).
strength_weakness: strength_weakness: - even though the authors claim the performance is "mixed" or better "in some cases", the proposed method performs significantly worse than most multi-step approaches across almost all datasets (see Table 4). In almost all other experiments (Table 3, 5, 6) the proposed method performs equal or worse than other multi-step approaches.
- I would then expect the proposed method to highlight a trade-off between performance and computational speed, but no thorough runtime comparison study is performed.
questions: questions: 1. as I understand the method, the main advantage comes from fast inference. Since the performance compared to other methods are worse or equal in most settings (except highly multi modal distributions), I believe it is important to show a runtime comparison with SOTA approaches. In particular, the very recent method in [1] is particularly flexible and much faster than other concurrent work. How does concurrent work compare with the proposed method at fixed performance (e.g. for a given log likelihood value) in terms of runtime (both training and inference)?
2. what are the numerical log likelihood values relative to the plots in Figure 6? how do they compare to other methods?
3. the proposed method relies on two fundamental approximations: estimation of the trace for the Jacobian determinant and approximate inverse of the learnt transformation. In the appendix the authors argue how to reduce the variance by rescaling v and that this helps particularly in low dimensions. I would be interested to see empirically how the variance scales with the dimensions in some experiments
4. given the several losses used in optimization (eq. 22) I was wondering how well is the inverse approximated in practice, which is the second fundamental approximation of the model. Also, the additional two losses kind of suggest that it is hard in practice to learn the approximate inference, which is however fundamental for the method to work. The authors also highlight that the reconstruction error plays a crucial role to learn the correct log likelihood. What are in practice the reconstruction error values (and the other loss terms) in the experiments you have run?

As a general advice, I would suggest the authors to further investigate how and why the proposed approach performs better when the learnt densities are multi modal. I think this is a very promising aspect - and from experience maybe the true limitation of bijective layers - but it is not explored in the paper

[1] Chen et al., Flow Matching on General Geometries, ICLR 2024

[2] Draxler et al., Free-form Flows: Make Any Architecture a Normalizing Flow, AISTATS 2024
limitations: limitations: Limitations are not at all clearly stated as such in the paper. Some limitations are referenced in the checklist but their description is vague and some references are missing (e.g. the checkerboard dataset. Without a comparison is also hard to grasp how much of a limitation this is). I would suggest the authors to include a paragraph concerning limitations, where limitations are mentioned and discussed.

Comments: title: Thank you for the rebuttal
comment: I would like to thank the authors for the detailed rebuttal.

As reviewers m8ZT and RvPv I was initially concerned about the reconstrcution error and I'm happy to see that in practice it is negligible. I now think about the proposed model as a fast alternative to competing methods and I think it is particularly relevant among single step methods. I still believe, however, that the comparison with multi step methods is not very convincing.

I will adjust my score accordingly.
Thank you for the response, I agree that small reconstruction errors might mitigate the likelihood issue and think that identifying how they are related is out of the scope of this paper.  Also thank you for reporting the W2 distance of your method on the various datasets, this is a great item for future works to compare against.  Since you have addressed all of my concerns, I'll update my score.
Thank you for the clarification and the continued interest! We do think that the reconstruction losses validate reporting likelihoods. Below, we provide more details why reconstruction errors are useful to detect the described behavior, as well as complimentary experimental evidence that no important regions are left out by M-FFF.

Let us consider the reviewer’s example of a uniform distribution, as this is an instructive case. Suppose that the model learns a uniform distribution, but leaves out a region $T$. Then we agree, it would result in a spurious likelihood. However, the reconstruction loss notices: A point in $T$ is moved by the encoder-decoder pair into $V \\setminus T$ (since the generator only produces points on $V \\setminus T$). Thus, all points in $T$ are moved into $V \\setminus T$, increasing reconstruction loss. The reconstruction losses consistently lower than 1 by three or more orders of magnitude indicate that no such left-out region exists or it is negligibly small.

Experimentally, we additionally verify that **every test data point has a generated point close by** on all $n=2$ manifolds (sphere, bunny, torus). To test this, we divide the manifold into bins and confirm that in each bin with a test data point, M-FFF also generates samples. We test this by sampling often enough from the trained model, until there is a sample in each bin. We choose natural bins for each manifold (icosphere, subdivision of original trimesh, rectangular grid of polar representation) and make the bins small (bin area $1/100,000$ compared to total manifold area). We confirm for every run for every $n=2$ experiment that in every bin populated by test data M-FFF also generates samples.

The problem with a sample-based metric is that competitors have not made their pretrained models available, so a thorough comparison would not be possible. We propose adding sample-based evaluations to the camera ready revision for future work to compare against. As a preliminary metric, we report the Wasserstein-2 distance for the $n=2$ datasets above:
| Dataset | Wasserstein-2 distance |
|---------------|---------------|
| Bunny (K=10)  | 0.09 ± 0.02   |
| Bunny (K=50)  | 0.046 ± 0.006 |
| Bunny (K=100) | 0.026 ± 0.007 |
| General       | 0.21 ± 0.04   |
| Glycine       | 0.32 ± 0.05   |
| Proline       | 0.51 ± 0.05   |
| PrePro        | 0.47 ± 0.04   |
| Flood         | 0.047 ± 0.010 |
| Fire          | 0.072 ± 0.027 |
| Volcano       | 0.249 ± 0.035 |
| Earthquake    | 0.068 ± 0.022 |
I appreciate the authors addressing my concerns.
- The authors provided the specific limitations on dimensions. 
- The authors were honest in their responses to my other questions.
- I strongly recommend that the authors include the content they provided in the updated version.

Overall, thank you for the efforts in addressing my concerns, and I have decided to adjust my recommendation to a score of 7.
Thank you for your rebuttal comments; in particular, the comments regarding when reconstruction error does not drop to a satisfactory level. I maintain my original score and believe this paper should be accepted.
Thank you for the response and for clearing up some of my concerns.  I appreciate adding reconstruction losses to the text however I am still concerned about using comparing likelihoods from the proposed model to those produced by exact flows.  My concern stems from the possibility that if the model does not map to the entire data space, then its values for probability density can be inflated compared to a model that maps to the entire data space.  For example, a uniform distribution over a space with volume $V$ has a density of $\frac{1}{V}$, but a uniform distribution over a subset of the region, with volume $V - T$ where $T>0$, has a larger density of $\frac{1}{V - T}$.  It is possible that something similar is happening in your experiments which is why I would like to see sample based evaluation.  Please let me know if you do not agree with this reasoning.