Score-Based Density Estimation from Pairwise Comparisons
Authors: Petrus Mikkola
Luigi Acerbi
Arto Klami
Venue: ICLR 2026 Poster (2026)
Keywords: score-based methods
pairwise comparisons
density estimation
elicitation
random utility models
tempering
Primary area: probabilistic methods (Bayesian methods, variational inference, sampling, UQ, etc.)
TL;DR: We show how to estimate densities solely from pairwise comparisons. We establish a relationship between the target density and tempered density of the preferred choices, and provide a score-based method that recovers the target density.

Abstract: We study density estimation from pairwise comparisons, motivated by expert knowledge elicitation and learning from human feedback. We relate the unobserved target density to a tempered winner density (marginal density of preferred choices), learning the winner's score via score-matching. This allows estimating the target by `de-tempering' the estimated winner density's score. We prove that the score vectors of the belief and the winner density are collinear, linked by a position-dependent tempering field. We give analytical formulas for this field and propose an estimator for it under the Bradley–Terry model. Using a diffusion model trained on tempered samples generated via score-scaled annealed Langevin dynamics, we can learn complex multivariate belief densities of simulated experts, from only hundreds to thousands of pairwise comparisons.

Published as a conference paper at ICLR 2026
SCORE-BASED DENSITY ESTIMATION FROM PAIRWISE COMPARISONS
Petrus Mikkola
Luigi Acerbi
Department of Computer Science University of Helsinki
Department Computer Science University of Helsinki
petrus petrus.mikkoon www.kkkkkkklllle.omammaa @iiaaehhelsssii................e
luigi luigi.acerbll www.baeeeeeeeiiieeeeeeeeeeeeeeee
Arto Klami'
Department of Computer Science University of Helsinki arto, k ami @he sinki. fi
ABSTRACT
We study density estimation from pairwise comparisons, motivated by expert knowledge elicitation and learning from human feedback. We relate the unobserved target density to tempered winner density (marginal density of preferred choices), learning the winner's score via score score-matching matching. This allows estimating the target by "de-tempering" the estimated winner density's score. We prove that the score vectors of the belief and the winner density are collinear, linked by a position- dependent tempering field. We give analytical formulas for this field and propose an estimator for it under the Bradley-Terry model. Using a diffusion model trained on tempered samples generated via score-scaled annealed Langevin dynamies, we can learn complex multivariate belief densities of simulated experts, from only hundreds to thousands of pairwise comparisons.
INTRODUCTION
Several complementary techniques, from flows (Rezende & Mohamed, 2015
Lipman et al. 2023) to diffusion models (Ho et al., 2020), can today efficiently learn complex densities p(x) from examples learn accurate densities even over high high-dimensional
p(x). With sufficiently large data, spaces, such as natural images (Rombach et al. 2022). While challenges persist in the most complex cases, these models have achieved a high level of performance, proving sufficient for many tasks.
We consider the fundamentally more challenging problem of learning the density not from direct observations but solely from comparisons two candidates. Given and x' that from the target p(x) but rather from a distinct sampling distribution X(x) satisfying suitable regularity conditions, the task is learn p(x) from triplets (x, x', > x'). The last entry indicates which alternative has higher density (the winner point). Being able to do this enables cognitively easy elicitation of subjective beliefs of an individual over random vectors. The canonical use-case is encoding expert knowledge into statistical models as prior information, with established literature in statistics dedicated to this problem of prior elicitation (O [O'Hagan] 2019: Mikkola et al.) 2023). Here the belief is typically over relatively low-dimensional space, but it needs to be inferred from very limited number of observations to keep the expert effort manageable. Recently, elicitation tools have been increasingly used to quantify large language model (LLM) knowledge in probabilistic terms (Capstick et al. 2025 Requeima et al. (0024), for instance, to evaluate calibration or to use them as probabilistic cognitive models (Binz & Schulz 2024) or forecasting models (Halawi al.) 2024). The current methods require dedicated techniques for direct prompting of probabilities, samples (Requeima al. 2024) or moments (Capstick al. 20225), whereas our formulation only requires comparative queries that LLMs can readily answer. Finally, the problem setup is also related to learning from human feedback (Ouyang et al. (2022), in particular to learning an implicit preference distribution for generative model (Dumoulin et al. (024)..
sampled
"qqual contribution.

---

Published as a conference paper at ICLR 2026
my belief
B X
"more probable"
(a) Problem setup
(b) MWD
(c) Tempering field
(d) 'Tempered' MWD
Figure I: (a) Problem setup. An expert holds a subjective belief over a parameter space, such as the likely hyperparameters of learning algorithm (e.g. learning rate and weight decay), and can answer questions like "Do you expect configuration A or B to work better?". We learn their belief as a density, to be used e.g. as a prior distribution for finding optimal hyperparameters. (b)-(d) Density estimation from 200 uniformly sampled pairwise comparisons, with the target density shown as heatmap. (b) Samples and the score field at an intermediate noise level a, for a diffusion model trained on the (winner, loser) pairs to model the marginal winner density (MWD). (c) Estimated tempering field. (d) Samples from the score- score-scaled annealed Langevin dynamics with the MWD score and a tempering field estimate. Samples align well with the target density, demonstrating the fundamental relationship between the scores of the estimable MWD and the latent target (belief density).
Recently, Mikkola et al. (2024) proposed the first solution for this problem, learning normalizing flows from pairwise comparisons and rankings. We propose an improved solution that also uses random utility models (RUMs: Train 2009) for modeling the preferential data and is inspired by their idea of relating the target density p(x) to tempered version of the distribution of winner points, Pi (x), for some tempering parameter T 2 1. Since we have samples from (x). this relationship leads to practical algorithms once is estimated. In contrast to their empirically motivated heuristic link, we characterize this connection in detail and provide an exact relationship between the scores of p(x) and (x). Since the relationship holds for the scores, it natural also switch to solving the problem with score-based models (Song& & Ermon 2019 Song et al.) (021). instead of flows This brings additional benefits, for instance in modeling multimodal targets, and we empirically demonstrate a substantial improvement in accuracy compared to Mikkola et al. (2024). While they could learn densities from modest number of rankings, they needed additional regularization to avoid escaping probability mass (Nicoli et al.) 2023). Moreover, their best accuracy required more informative multiple-item rankings. In contrast, we focus solely on pairwise comparisons, which are easier to answer and more reliable (Kendall & Babington Smith 1940, Shhh & Oppenheimer 2008). and widely used in Al alignment (Ouyang et al. [0022] Wallace et al. 2024).
Denote by Px-xx (x, x') the joint density of the available data, encoding the preferred candidate in the order of the arguments. The marginal winner density (MWD), denoted by (x). is obtained where A(x) is the
as its marginal sampling density of the (independent) candidates. Our main theoretical contribution is a novel, exact relationship between the target p(x) and the MWD Par (x) in terms of their scores: up reparameterization of the space, we have V loggpf p(x) T(x) V log (x). Critically, 7(x) is not constant but position-dependent tempering field. This implies we can perfectly recover p(x) from the estimable (x) with score- score-based methods if the tempering field is known. We prove this foundational relationship for the popular Bradley-Terry model (Bradley & Terry 1952, Touvron et al. 2023) and an exponential noise RUM, providing explicit formulas for the tempering fields.
(x) = Px>xx (x,x'xxxxx))) x')dx' x [P(x x')x(x)Xxx')
Our second contribution a practical algorithm derived from our theoretical insights. First, we pro- pose to model the preference relationships by estimating the score of the joint density Px- (x.x). then train a continuous-time diffusion model (Karras et al. 2022) to recover the MWD by marginal- izing it. Building on the ideal tempering field under the Bradley-Terry model, we estimate the tempering field 5 rxx) by using the analytical formula with importance samples from the trained MWD model and simple density ratio model trained on the pairwise comparison data. Finally, we sample from the belief density p(x) by running score-scaled annealed Langevin dynamics (Song & Ermon 2019) with the MWD score and rxx). Fig. [ illusttrates our approach.
2

---

Published as a conference paper at ICLR 2026
2 BACKGROUND
2.1 DENOISING SCORE MATCHING AND ANNEALED LANGEVIN DYNAMICS
The (Stein) score of probability density function p(x), denoted Vx log p(x), is vector field pointing in the direction of maximum log-density increase. Score- base- generative methods approximate this score. They typically start by defining family of perturbed densities Pa (x) by convolving p(x) with noise at varying levels a > O; for example, Px( = p(x) (xx 0, oT), where denotes convolution. A neural network So 5x(x,77 07) with parameters 0 is then trained model the score of these perturbed densities, Vx log Pa (x). This score network Sa is commonly trained through demoising score matching (Vincent 2011. by minimizing the objective:
(())
Bxxp, (o)) V (x sxx
(1)
Here, & is noisy version of a clean sample x, generated via the perturbation kernel po (x|x) (e.g. an isotropic Gaussian (8; x, I²)).. The network is trained predict the score of Per by minimizing the objective in Eq. where the perturbation kernel is typically tractable. The function {((r provides a positive weighting for different noise levels. The noise levels g are drawn from distribution Purin (o)
following either discrete, often uniform schedule (or) (Song & Ermon 2019), or a continuous one Karras al. 2022).
Once trained, (x a)) enables sampling from an approximation of p(x). One prominent method, besides reverse diffusion processes (discussed later), annealed Langevin dynamics (ALD) (Song
& Ermon 2019). ALD starts with samples XT from broad prior (e.g., (<( 0, o I)) and iteratively refines them. runs L steps of Langevin MCMC per noise level o, along a decreasing schedule max
01 = mmin:
x (1-1) + (x (1-1) .x) + 26, n;
1=1,2,...,
(2)
with step size > 0 and n! ~ ~(O,II, I). Fort < T,
x!" (0)
Xi41 Under ideal conditions (L oC.
(L)
(L) approximates sample from Pa (x) ~ p(x) (Welling & Teh 2011
O, accurate Sq), xi
2.2 DIFFUSION MODELS
continuous-time diffusion model describes a forward process that gradually transforms a data distribution p(x) into simple, known prior distribution (e.g., a Gaussian). This process is often defined by forward-time stochastic differential equation (SDE) (Song et al. 202I
dx = f(x, l)dt + g())db,
where b is Brownian motion, f(x, t) is the drift coefficient, and g(t) is the diffuuson coefficient. If x(0) ~ p(x) (the target density), its time-evolved density is (x)).. If f is an affine transformation, then the transition kernel p(x(f) x (0) is Gaussian and for a sufficiently large T > 0, the marginal distribution Pr (x(()) becomes pure Gaussian, such as No I) or ( II).
The forward process can be reversed to generate data. Starting from a sample xy pr (T(x). one can obtain a sample Xy ~ p(x) by solving the corresponding reverse- time SDE (Anderson 1982);
dx (J(x, 1) (1) log (x)) di + g()ddb,
where b is Brownian motion with time flowing backward from T to 0. Alternatively, samples can be generated by solving the deterministic probability flow ODE (Song et al. 2021).
dx
f(x - (1)
dt.
Both reverse methods require the score function Vx log Pt (x), typically approximated by a trained score network, sx t) or So Sx .r) if parameterized by noise level o.
The Elucidating Diffusion Models (EDM) framework (Karras et al.
2022 2024a) parametrizes the diffusion process directly using the noise level o [miin mmax rather than abstract time
t. This can be achieved by assuming gt) /2t and f(x, t) = 0, and using the initial condition Xy AO, 2 I) for some fixed, sufficiently large Fmax > 0. The perturbed density can be written (x) Po (x) p(x) A Nxxx; I
The score network sy (x. o)) is trained via denoising score matching (Eq. []. Sampling is done by solving the stochastic reverse diffusion SDE (Eq. 3 the deterministic probability flow ODE (Eq. 4).
3

---

Published as a conference paper at ICLR 2026
2.3 RANDOM UTILITY MODELS AND DENSITY ESTIMATION FROM CHOICE DATA
In the context of decision theory, the random utilityy model (RUM) represents the decision maker's stochastic utility U as the sum of deterministic utility and a stochastic perturbation (Tain 2009).
(xx) u(x) + (xx),
where u X' -> R is a deterministic utility function, W is a stochastic noise process, and the choice space X' is compact subset of Rdd Given a set C C X of possible alternatives, the decision maker selects an alternative x E C. by solving the noisy utility maximization problem: x arg ea U (x'). Pairwise comparison is the most common form of choice data and corresponds to assuming that the choice set contains only two alternatives, C= {x, x'}. The decision maker chooses x from C, denoted by x', if u(x) + w(x) > u(x')+ + w(x') for a given realization w of W. It is often assumed that W is independent across x, leading to so-called Fechnerian models (Becker et al., 1963), where the choice distribution conditional on C reduces to F(u(x) =(x)),
with F' denoting the cumulative distribution function of W(x) (x)..
Psychophysical experiments suggest that human perception of numerical magnitude follows a log- arithmic scale (Dehaene) 2003 Assuming a RUM with utility function u(x) log p(x), the model's noise becomes additive in the log-transe log-transformed beliefs. In this paper, we consider two RUMs, explicitly including the noise level, as it is crucial for identifying p(x First, we study the generalized Bradley-Terny model (Bradley & Terry 1952) with W ~ Gumbel s), which induces the conditional choice distribution l'Logisti((O.e)) (x() u(x')). Second, we consider the exponential RUM with W ~ Exp(s), which yields a heavier-tailed conditional choice distribution Lapplace(0,11 s) (w(x) u(x))).
Under these assumptions, the density estimation task is instance of expert knowledge elicitation OOHaaga 2019: Mikkola et al.] 2023). and it is closely related to (probabilistic) reward modeling Leike al. [2018] Dumoulin et al. (024). Expert knowledge elicitation infers expert's belief as
probability density p(x) using only queries they can reliably answer, such as requests for specific quantiles or statisties of p(x) (O'Hagan) 2019; [Bockingg et al. 2025) or comparisons like here. Recently, [Dmmoulin et al. (2024) reinterpreted reward modeling by referring to the target distribution as the "implicit preference distribution" and treating the reward as a probability distribution to be
modeled.
3 BELIEF DENSITY AS A TEMPERED MARGINAL WINNER DENSITY
Let p(x) be the expert's belief density. We assume the expert's is choices follow RUM with utility function u(x) = log p(x). They observe two points independently drawn from the sampling density A(x), and the expert chooses one of the points. We denote the probability density of that point by Pw (x) and refer to it as the marginal winner density (MWD). By marginalizing out the unobserved loser x' in a pairwise comparison where x is preferred (x > x'), (x) can be expressed as 2X(x) f F(log((x) - log ((xxx)))x''ddd )dx (Mikkola et al. 2024).
While Mikkola et al. (2024) empirically showed that p(x) resembles a tempered version of p.r ((x...(() R [x]] for some constant r), this relationship was not formally analyzed besides the theoretical limiting case of selecting the winner from infinitely many alternatives. In this section, we establish more fundamental connection. We demonstrate that under two RUMs-the Bradley- Terry model and the exponential RUM-it is possible to find a tempering field 7(x) such that V Plogp(x) 7(x) V log Pw (x) up to reparameterization of the space. This key relationship implies that, in principle, p(x) can be precisely recovered from (x) using score-based methods if 7xx) is known. This finding motivates leveraging acore-macching techniques for estimating the belief density. To analyze such score-based relationships and evaluate approximations, we use the Fisher divergence, which quantifies the difference between two distributions based on their scores:
F(p.g) = [VVlog(x) V log g(x) p(x)dx.
EExact identification of p(x) requires knowing both the correct noise family and noise level. We consider the noise he known; see Appendix F.3 for validation of robustness this choice.
4

---

Published as a conference paper at ICLR 2026
3.1
TEMPERING FIELD
Consider a tempered probability density p(x) constructed from another density (xx) using a tempering constant 0:
9" (x)
p(x) Jx q" ('('')dx''
For such densities, the relationship between their scores is given by the product rule as
V plg)(()) 7V Vologg (x)) log (())VT rV log (xx).
The score of the tempered density becomes directly proportional to that of the original density, i.c., the two score vectors are collinear. Inspired by this relationship, we define a more general concept. We call function X' (0, 20) a tempering field between and q if their scores satisfy the following relation almost everywhere for x E X::
Vlog = T Tx)) log (()).
(5)
This implies the scores collinear, with rxx) position-dependent scaling. The tempering field thus describes a localized, score-level tempering.
3.2 TEMPERING FIELDS UNDER RUMS
With our theoretical framework in place, we analyze the relationship between the belief density p and the MWD Pw in terms of tempering fields for RUM models with utility logp. We prove our main results for both the Bradley-Tert model and the exponential RUM, with W ~ Gumbel(0, s) and
W Exp(s). The treatment of the latter RUM is deferred to Appendix A
To facilitate theoretical analysis, with no loss of generality, we assume uniform sampling distribution A throughout this section, to remove the tilting of MWD) Pwo (x). We then address non-uniform À by reparameterizing the space so that it becomes uniform on a hypercube. The invariance of the RUM under the space reparameterizing is derived in Appendix [ The diffusion model is trained in the transformed space, and the generated samples are mapped back to the original space with the inverse transformation. We use the Rosenblatt transformation, which requires the conditional distribution functions of X (Rosenblatt] 1952). here assumed to be either known (e.g., when À is Gaussian) or numerically approximated. Other transformations, e.g. ones based on copulas or normalizing flows trained on samples from À (i.e., the combined data of winners and losers), could be used as well.
Under the following assumptions, a tempering field exists between the belief density and the MWD: Assumption 1. supp(p) C supp(A).
Assumption 2. À is a uniform density over X.
Assumption 3. is smooth, with Vp + O almost everywhere.
Theorem 3.1. Assume W ~ Gumbel(O, s). A tempering p and the MWD Prs and given by the formula,
field field (xx) exists between the belief density
Jx I+r, (x,] dx' To (x(xxx)
T(x)
(6)
Sx (I+rx(xx]]> dxx
where r,(xx x')= = p² (x')p - (xx
(x) is the teeel//sseeepe 1/8-tempered density ratio.
Proof. Our constructive proof derives a scalar field that satisfies the defining relation. Specifi- cally, for any fixed x E X, direct manipulations yield a scalar 7(x) >> 0 such that Vx logp(x)
T T(x) Vx Pw. = 0. See Appendix B for the full proof.
Fig. 2 illustrates the tempering field in Theorem 3.1 using s = /6/m² (unit variance noise). There exists specific invariance relationship between the tempering and the noise scale. Specifically,
Tp.s denotes the tempering field of RUM under the belief density p and noise level 8 > 0, then by the tempering field theorems it is clear that for any exponent > 0: Tp _ and Tpm = ST 2 11' where the tempering fields are of the exponential RUM and the Bradley-Terry model, respectively.
5

---

Published as a conference paper at ICLR 2026
(a) Score fields
(b) Tempering fields (estimate and ground-truth)
Figure 2: Illustration the relationship logp(x) T(x) V log Pw when p is Twomoons2D (Stimper et al. 2022) and is uniform. (a) The score of (red arrows) and the score of pw (orange arrows) under the Bradley-Terry model, scaled for better visualization. (b) The estimated tempering field (xx from 200 pairwise comparisons (left, Section 4.2) and the ground- -tuth (right, Theorem 3.1 Due to the collinearity of the scores, the red arrows equal the pointwise product of the orange arrows and the tempering field, which can be estimated (with underestimation in this example).
3.3 ON THE CONSTANT TEMPERING APPROXIMATION
Even though our method directly estimates the full tempering field, our theory also sheds light on methods assuming constant tempering, such as Mikkola al. (2024). It allows us to establish three quantities of interest related approximating with a constant-tempered version of q: (i) the optimal constant tempering T* >> O which minimizes F(p, q²') (ii) the approximation error F(p, q²) for any constant > 0, and (iii) the approximation error F(p, q" for the optimal constant tempering. Proposition 3.2. Assume that there exists a tempering field r(x) between and The optimal tempering constant T" arg min,>0 F(p, 4") can written as,
Ex-p (x(X)
where the stochastic weight w 2 0 is given by w(X IV log x(X) Ex-p log (())²)
Proof. By the Leibniz integral rule,
F(p.a") = / 2 (1lIV (V log (x),, V Vlogg(x)))) p(x)dx. 0r
The divergence is quadratic in t and the critical point is the global minimum, and by assumption (V log (xx), V log ((x)) r(x) V log ( Algebraic manipulation gives the result.
The approximation errors can be quantified in terms the tempering field. Proposition 3.3. Let r(x) be a tempering field. For any 0 it holds that
F(p.4") = Ex-p (Ir (I--r(X)) V (
Further, when r* is the optimal tempering, we have
F(p.q"') Exap (1) (Exap (7(X) V/ log g( x(x)11))
Exap (I log 1
Proof. See Appendix

---

Published as a conference paper at ICLR 2026
4 SCORE-BASED DENSITY ESTIMATION FROM PAIRWISE COMPARISONS
Building on Section 3 we now introduce our score-based density estimator for eliciting the belief density from pairwise comparisons. The method has two components. First, train a diffusion model on the joint distribution of winners and losers using a masking scheme that ensures its marginal, that is MWD Pwu (x), can also be evaluated. Second, under the Bradley-Terry model, we provide a simple procedure to estimate the tempering field r(x) and use it to de-temper the score-based estimate of the MWD. Details of both steps are explained next. The sampling distribution A A(x) is assumed known, and we reparameterize the space to make it uniform, as explained in Section 3.2
4.1 MODELING THE MWD
Our goal is to learn the perturbed score model of the MWD V log(p=((x) * <xx; 0, D) We want to utilize all samples, both winners and losers. To do so we simultaneously learn the marginal Pw (x) = J Px--x (x, x') 2')''', and the full joint Px>xx (x, x') from the concatenated data of winners and losers. To learn the marginal, during training, half of the time we randomly mask x' and consider the score only with respect to
We parametrize the score model as 50(x, x', ,. joijn,,, temp). where joint E {0,, 1 and temp E {0, 1) conditioning flags: (a) When joint = 1, the network takes both and x' as input and models the score of the joint distribution, V loog Px-xx (x, x'). (b) When joint = 0, the loser x' is masked (replaced with noise), and the network models the MWD score V log Pie Prc (x). The flag temp then controls whether the output is scaled by the tempering field: setting temp = yields an approximation to 7xx) log Pw (x) V log p(x).
This parametrization allows us to train a single score network on both winners and losers via denoising score matching (Eq. D while still enabling sampling from the belief density. During training, we randomly mask the loser with probability 0.5 (see Appendix C.T for details). The MWD could also be estimated directly from winners alone, but this would ignore that losers carry information about where winners are less likely to be. We demonstrate the value of modeling the joint distribution empirically in Appendix while also confirming that we can accurately marginalize the joint model.
We stay as close as possible to the EDM-style diffusion model (Karras et al.) 2024a). We use the perturbation kernel (x x) Nixx; x, I²I), which aligns with EDM and defines forward diffusion process from Jmin to mmax. where Pamii (x) Ri p(x) and Po (x) 2 NO o diax A
detailed description is provided in Appendix C.2
Algorithm summarizes the method.
4.2 TEMPERING FIELD ESTIMATION
The tempering field under the Bradley-Tenry model (Bq. 6) has particularly convenient form as it depends only on the belief density ratio xxx := p(x')/p(xx and the RUM noise level > 0. Note that this ratio is different from what the phrase density ratio often refers to; this is the ratio of the same density for two inputs, not a ratio of two densities for the same x. It does not depend on the normalizing constant of the belief density and is hence straightforward to estimate via estimation (MLE), assuming careful regularization. We train a simple estimator for r(x, x') by maximizing the Bradley-Terry model likelihood of the pairwise comparison data. If we parametrize the density ratio (or its logarithm) as a neural network re, the parameters @ can be optimized by minimizing the loss L(Ø) Softplus(log (x x))/s,, where x and x' are winner and loser points, respectively. The Softplus loss comes from the assumption of W ~ Gumbel s).
We plug the learned ro into the integrals in Eq. 6 where the integrals are computed using importance sampling with the MWD model acting as the importance sampler. The resulting plug-in Monte Carlo estimator of the tempering field is consistent but biased. Similar biased ratio estimators have been used in self-normalized importance sampling (Owen) 2013) and in every-visit off-policy value estimation in reinforcement learning (Sutton et al., 1998) due to favorable variance properties. See Appendix C..5 for more discussion the estimator. For details evaluating the importance weights, which correspond to the (reciprocal of) density of the MWD diffusion model, see Appendix C.7 Algorithm summarizes the estimation procedure.
7

---

Published as a conference paper at ICLR 2026
4.3 BELIEF DENSITY SAMPLING
Given the perturbed MWD score network (x(x, 2 V log[pm((x) * (x 0, o²1)] and the estimate of the tempering field 7(x), we can draw approximate samples from the belief density p(x) using the score-scaled ALD. Specifically, we iteratively run Eq. with the score T T( (x))ssx(xo)))))) o). The tempering field relation (Eq. 5 shows that for = 0 this procedure would be equivalent to sampling from p(x) and ALD is theoretically valid at the small-noise limit (Welling & Teh, 2011}, making it an appealing choice as sampling algorithm. However, for > the algorithm not exact. We show that it works empirically well (Section 5). but characterization of the approximation error remains as future work.
Algorithm 1 Full algorithm
Algorithm rxx)
require:
choice data D = {[x, x X; > x))==1
require: noise level s, D, 8g initialize:
output: samples from the belief density or a trained diffusion
model for it it
train re x') R3 p(x')
p(x) as
train (x, x', a, joint) using DSM (Eq. ! D
MLE of 'D (Softplus loss)
50% prob.: set joint = and mask x' with W
sample m points X with den-
50% prob.: set noise schedule (or)
sities pw (X) using
initialize 7(x) given D and Sp
input:
sample D using -sx)-e ALD with the score
r (ro(x, X))
So x, x", (ot) optional
return:
147, (X,)
train SBuwp (x, G, temp) using DSM on D
(1+r)) px(
return: DD* or S0wen (x, o, temp=1)
5 EXPERIMENTS
We evaluate the method on synthetic data generated from a RUM. We then consider an experiment where a large language model (LLM) serves as a proxy for human expert, demonstrating the method's applicability in settings where data does not follow a RUM model. Our experimental setup closely follows that Mikkola al. (2024), with the key distinction that we only query pairwise comparisons, not considering the easier case of ranking multiple candidates. We empirically compare against their fowed model, using their implementation, as the only previous method for the task. We consider two variants of our score-based method: score-7(x), which uses the full tempering field (Section 4. and score-r", which uses a constant tempering estimated via Proposition 3.2 This allows direct quantification of the importance of modeling the whole field.
Setup and evaluation For a d-dimensional target, we query 1000d pairwise comparisons to ensure reliable comparison between the methods but remaining well below the surge-sample regime typical for the diffusion model literature. For d < 4, the sampling distribution À is uniform. For d > 4, is a diagonal Gaussian (Gaussian mixture in Mixturegaussians 10D) 10D) centered at the target mean, with a variance three times that of the target's The simulated expert follows the Bradley-Terry model with utility given by log p. where the belief density p varies in dimensionality, modality, and detail. We set the noise level 8 / /6/7² (unit variance). As the diffusion model, we adopt an EDM-style architecture with MLP score network, implemented top of Karras al. (2024b). For further details, see Appendix C and We assess performance qualitatively by visually comparing 2D and 1 12 marginals of the target density and the estimate, and quantitatively using two metries: the Wasserstein distance and the mean marginal total variation distance (MMTV; Acerbi 2020). Results reported as means and standard deviations over replicate runs.
Experiment 1: Synthetic low dimensional targets with uniform sampling We consider low- dimensional synthetic target densities that may exhibit non non-trivial trivial geometry or multimodality. The set (top) shows that the score-based method is clearly superior for all targets, with at least 50% (Wasserstein) and 25% (MMTV) reduction of error
includes five target distributions, see Appendix E. Table |
For non-uniform sampling distributions, we transform the space with Rosenblatt transform such that the sampling pdf becomes uniform, with an appropriate Jacobian transformation to the densities.
8

---

Published as a conference paper at ICLR 2026
Table 1: Density estimation from pairwise comparisons: score-7* and score-r(x) denote our
methods with constant and varying tempering fields, respectively, and flow refers to (Mikkola et al. 2024). Bold indicates the best method, and underline indicates results that are not significantly worse (paired two-sided Wilcoxon signed-rank test, p > 0.05).
wasserstein (L)
MMTV (4)
p(x)
flow
score-r*
(xx
flow
score-T"
score (xx
Onemoon2D
1.37 (=0.03)
0.44 (#0.13)
0.37 (#0.14) 0.54 (=00000
0.25 (#0.06)
0.22 (#:0.06)
Twomoons2D
1.29 (=0.06)
0.54 (#0.14)
0.44 (#0.09)
0.53 (#0.01)
0.15 (#0.03)
0.14 (#0002)
Ring2D
0.87 (=0.03)
0.40 (#0007))
0.39 (=0.08)
0.40 (#0.01)
0.27 (:0001))
0.26 (#:0.01)
Gaussian4D
6.12 (+0003)
1.69 (4-0.22)
1.40 (#0.24)
0.72 (+0.01)
0.41 (+0005)
0.44 (++0.08)
Mixturegaussians4D
3.75 (#0.02))
1.23 (4-0.66)
1.09 (#0.09)
0.53 (+0.011
0.27 (+00011)
0.22 (++0002)
Stargaussian6D
2.25 (=0.02)
T.55 (#0.04)
1.28 04)
0.79 (20.00)
0.T8 (20:02)
0.16 (#0.01)
Mixturegausiaas
1.41 (#0.01)
1.10 (#00122)
1.33 (#0.11)
0.19 (#0.00)
0.14 (#0.02)
0.26 (#0.00)
Gaussian16D
5.50 (+0.03)
5.00 (#0.07)
5.00 (#0.08)
0.16 (+0.00)
0.13 (+0.00)
0.13 (+0.00)
compared to the flow method. most cases, using the full tempering field Txx) is better than relying on the best constant tempering, but we outperform the flow-based method even when restricted to constant tempering as in their approach. Visual inspection (Fig. 3 (a-b) and Figs. F.1 F33 confirms this. The flow method captures the density relatively well but clearly overestimates the low-density regions, whereas our estimate is essentially perfect here.
Synthetic targets with Gaussian sampling For higher-dimensionta experiments we replace uniform sampling with more concentrated Gaussian sampling, since otherwise the proba- bility that both and x' fall in low-density regions increases dramatically, making it impossible to learn p(x) well. The set includes three target distributions, see Appendix E.1 Table (bottom) shows that we again outperform the flow-based method in Wasserstein distance but in terms of MMTV the methods are closer. Visually, the score- score-based method usually gives sharper and better estimates (e.g., compare Fig. and
Experiment 2:
and F.7). . but suffers in terms of the MMTV metric due to occasional too-tight marginals resulting from overestimating the tempering field (e.g., Fig. F.8).
Experiment 3: LLM as a proxy for the expert To illustrate the method in real belief density estimation task without user experiments, we replicate the LLM experiment of Mikkola et al. (2024). except that we query 220 pairwise comparisons instead of 5-wise rankings. Using the data from (Mikkola et al. 2024), we uniformly sample 220 pairwise comparisons across all eight features and prompt the LLM for belief judgments. The LLM acts as a proxy for a human expert, providing its belief about what houses in California are like, restricted the features available in the California housing dataset (Pace& & Barry 1997). The LLM's s belief is inferred solely from the pairwise queries, without providing the LLM any direct access the data itself.
This experiment leverages the finding that L.Ms learn statistical features from the vast training data and can be queried about it (Brown et al.) 2020 Requeima et al. 2024). While the aim is infer the belief density of the LLM-which is unknown- our main goal here is to validate our method, rather than to analyze the beliefs of this particular LLM. We do this by comparing the belief estimate with the empirical data distribution: Similarity between the two suggests the elicitation method yields a reasonable belief estimate, and differences can be interpreted as possible biases the LLM might have. Fig. (c) shows clear distributional similarities; e.g. the marginals of AveRooms and MedInc exhibit similar shapes. See Appendix F22 for complete results, with Table F. quantifying the accuracy.
DISCUSSION
We proved the theoretical connection for two common RUMs but we believe it extends to other RUMs as well, although a closed-form tempering field is not guaranteed e.g. for the Mostelee model (Thurstone 1927) where the choice probability requires integration.
The difficulty of the belief estimation problem, naturally, depends heavily on the sampling distribution X(x). For example, when the support of p(x) is much smaller than that of X(x). it becomes
Foor a direct comparison with Mikkola 2024 we also use Claude Haiku by Anthropic 2024

---

Published as a conference paper at ICLR 2026
Leamed LLM
Colifornia housing dataset distribution
(a) Score-hased
(b) Flow
(c) LLM experiment
Figure 3: (a-b) Samples from score-based and flow estimates of Ring2D, with contours indicating the true density. Ring2I) illustrates an extreme where the score-based method clearly outperforms the flow method: the flow model oversamples the center of the ring, where the MWD also has moderate density, whereas the score-based method can downweight using the tempering field. (c) Cross-plot of the first three variables in the LLM expert elicitation experiment. Full cross-plot and comparison the flow method are shown in Figs. F. and F.11 The score-based method tends to generate Gaussian-lic like marginals in this extremely limited data setting.
nearly impossible to obtain sharp estimates of p(x), and the problem is further exacerbated in high- dimensional spaces. For this reason, we see potential in active learning methods that concentrate sampling in high-density regions of p(x). In this work, we assume that X(x) is given as it would be in many applications (c.g. the clicitation protocol in prior clicitation), but for active learning or learning beliefs from public preference data, additional density estimation step is required to learn X(x).
We demonstrated that low-dimensional targets can be learned from few hundreds of pairwise com- parisons (Fig. [ However, in extremely limited data regimes, say below 100d pairwise comparisons, the robustness of our method is not guaranteed without carefully tuning hyperparameters such as those of the diffusion model-aa class of models that notoriously difficult train in stable manner KKrrra er al. 2024bb Appendix B.5). Stability could likely be improved by adopting best practices from the field (Karras al. 2022) and incorporating recent advances in learning score models from limited data (Li al. 2022). Similarly, the tempering field estimate (Section 4.2) depends on the density ratio estimate To. which is sensitive the network's 's regularization. Misspecified values of the lg regularization for the network weights 0 or the RUM noise level s will lead to under- or overestimation of the tempering field, due to MLE struggling to determine the global scale and tails.
While our method is theoretically grounded and clearly outperforms the flow-based method of Mikkola et al.) (2024) in estimation accuracy, has some limitations. The flow allows faster sampling, only requiring a single forward pass, and also more efficient and stable evaluation of the density. Our method requires numerically solving the probability- flow ODE for evaluating the density, and is still open whether diffusion models provide reliable pointwise density estimates (Zheng et al. 2023).
The connection between the target density p and the MWD Pw may have applications beyond expert knowledge elicitation, such fine-tuning generative models (Wallace al. 2024) using pairwise data, perspective explored by Dumoulin et al. (024).. Specifically, when A(x) represents a pretrained generative model conditioned on prompt c, our theory suggests that training the MWD and tempering field on individual-level data (rather than pooled data) yields a probabilistic reward model reward(c, x) p(x[c). This conditional tempered MWD captures the distribution of samples (c.g., images) aligned with prompt from that specific individual's 's perspective.
7 CONCLUSIONS
Despite two decades active research on learning from preference data, following the pioneering works of Chu & Ghahramani (2005); Fürnkranz & Hüllermeier (2010), the question of how to learn flexible densities from such data has remained elusive. We established the missing theoretical basis by showing how the score of the target density relates quantities that can be estimated, enabling the use of powerful modern density estimators for this task. Our approach demonstrates superior performance recent flow-based solutions in representing human beliefs.
10

---

Published as a conference paper at ICLR 2026
ACKNOWLEDGMENTS
The authors acknowledge the research environment provided by ELLIS Institute Finland. The work was supported by the Research Council of Finland Flagship programme: Finmish Center for Artificial Intelligence FCAI, and additionally by the grants 363317 and 358980. The authors acknowledge support from CSC IT Center for Science, Finland, for computational resources.
Reproducibility statement The source code for reproducing the experiments is available at Section describes the method and its im- plementation, while Appendix C details the specific components and training settings required to replicate the experiments. Appendix provides the full proofs of the theoretical results presented in the main text.
http://gitubb.com/petuuuu http://hw....c.mmomommmmmmmoo//a
REFERENCES
Luigi Acerbi. Variational Bayesian Monte Carlo with noisy likelihoods. In Advances in Neural Information Processing Systems, volume 33, pp. 8211-8222. 2020.
Brian Anderson. Reverse-time diffusion equation models. Stochastic Processes and their Applications, 12(3):313-326, 1982.
Anthropic. The Claude Model Family: Opus, Somnet, Haiku, http:::/ http://wwwww.c..tton.....o..... Com/ ott.::/0.....a..a..........0...3
Card Claude
pdf 2024.
Model card.
Gordon M Becker, Morris H DeGroot, and Jacob Marschak. Stochastic models of choice behavior. Behavioral science, 8(1):41-55, 1963.
M Binz and E Schulz. Turning large language models into cognitive models. In Twelfth International Conference on Learning Representations, 2024.
F. Bockting, S. T. Radev, and P. C. Bürkner. Expert-elicitation method for non-parametric joint priors using normalizing flows. Statistics and Computing, 2025.
Ralph Allan Bradley and Milton Terry. Rank analysis of incomplete block designs: I. the method of paired comparisons. Biometrika, 39(3/4):324 345- 1952.
Tom Brown, Benjamin Mann, Nick Ryder, Melanie Subbiah, Jared D Kaplan, Prafulla Dhariwal, Arvind Neelakantan, Pranav Shyam, Girish Sastry, Amanda Askell, Sandhini Agarwal, Ariel Herbert-Vons Voss, Gretchen Krueger, Tom Henighan, Rewon Child, Aditya Ramesh, Daniel Ziegler, Jeffrey Wu, Clemens Winter, Chris Hesse, Mark Chen, Eric Sigler, Mateusz Litwin, Scott Gray, Benjamin Chess, Jack Clark, Christopher Berner, Sam McCandlish, Alec Radford, Ilya Sutskever, and Dario Amodei. Language models are few-shot learners. In Advances in Neural Information Processing Systems, volume 33. pp. 1877-1901, 2020.
Alexander Capstick, Rahul Krishnan, and Payam Barnaghi. AutoElicit: Using large language models for expert prior elicitation in predictive modelling. In Aarti Singh, Maryam Fazel, Daniel Hsu. Simon Lacoste-Julies Felix Berkenkamp, Tegan Maharaj, Kiri Wagstaff, and Jerry Zhu (eds.), Proceedings of the 42nd International Conference on Machine Learning, volume 267, pp. 6746-6777, 2025.
Ricky TQ Chen. torchdiffeq, 2018. URL http://pp///gttbb...c/////t//ttp
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud. Neural ordinary differential equations. In Advances in Neural Information Processing Systems, volume 31, 2018.
Wei Chu and Zoubin Ghahramani. Preference learning with Gaussian processes. In Proceedings of the 22nd International Conference Machine learning, pp. 137-144, 2005.
Stanislas Dehaene. The neural basis of the Weber-Fechner law: logarithmic mental number line. Trendy in cognitive sciences, 7(4):145-144 2003.

---

Published as a conference paper at ICLR 2026
Vincent Dumoulin, Daniel D. Johnson, Pablo Samuel Castro, Hugo Larochelle, and Yann Dauphin. A density estimation perspective on learning from pairwise human preferences. Transactions on Machine Learning Research, 2024. ISSN 2835-8856.
Harley Flanders. Differentiation under the integral sign. The American Mathematical Monthly, 80(6): 615-627, 1973.
Johannes Fürnkranz and Eyke Hüllermeier. Preference learning and ranking by pairwise comparison. In Preference learning, pp. 65-82. Springer, 2010.
Manuel Gloeckler, Michael Deistler, Christian Dietrich Weilbach, Frank Wood, and Jakob Macke. All-in-one iiaullationbascle inference. In Proceedings of the 41th International Conference Machine Learning, pp. 15735-15766, 2024.
Will Grathwohl, Ricky TQ Chen, Jesse Bettencourt, Ilya Sutskever, and David Duvenaud. Ffjord: Free-form continuous dynamics for scalable reversible generative models. In International Confer- ence on Learning Representations, 2019.
Danny Halawi, Fred Zhang, Chen Yueh-Han, and Jacob Steinhardt. Approaching human-level forecasting with language models. In Advances in Neural Information Processing Systems, volume 37, pp. 10426-50466) 2024.
Dan Hendrycks and Kevin Gimpel. 1606.0866.111 2016.
Gaussian error linear units (gelus). arXiv preprint
Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In Advances in Neural Information Processing Systems, volume 33, pp. 6840-6851, 2020.
Tero Karras, Miika Aittala, Timo Aila, and Samuli Laine. Elucidating the design space of diffusion- based generative models. In Advances in Neural Information Processing Systems, volume 35, pp. 26565-26577, 2022.
Tero Karras, Miika Aittala, Tuomas Kynkäänniemi, Jaakko Lehtinen, Timo Aila, and Samuli Laine. Guiding diffusion model with bad version of itself. In Advances in Neural Information Processing Systems, volume 37, pp. 52996-53021, 2024a.
Tero Karras, Miika Aittala, Jaakko Lehtinen, Janne Hellsten, Timo Aila, and Samuli Laine. Analyzing and improving the training dynamics of diffusion models. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pp. 24174-24184 2024b.
M. G. Kendall and B. Babington Smith. On the method of paired comparisons. Biometrika, 31(4/4): 324-345, March 1940.
Diederik P. Kingma and Jimmy Ba. Adam: method for stochastic optimization. In International Conference on Learning Representations, 2015.
Jan Leike, David Krueger, Tom Everitt, Miljan Martic, Vishal Maini, and Shane Legg. Scalable agent alignment via reward modeling: research direction. arXiv preprint 181.::111..07871 2018
Xiang Li. Yixiang Dai, and Qing Qu. Understanding generalizability of diffusion models requires rethinking the hidden Gaussian structure. In Advances in Neural Information Processing Systems, volume 37, pp. 57499-57533 2024.
Yaron Lipman, Ricky T. Q. Chen, Heli Ben-Hamu, Maximilian Nickel, and Matthew Le. Flow matching for generative modeling. In International Conference on Learning Representations. 2023.
llya Loshchilov and Frank Hutter. Decoupled weight decay regularization. In International Confer- ence on Learning Representations, 2019.
Petrus Mikkola, Osvaldo A. Martin, Suyog Chandramouli, Marcelo Hartmann, Oriol Abril Pla, Owen Thomas, Henri Pesonen, Jukka Corander, Aki Vehtari, Samuel Kaski, Paul-Christian Bürkner, and Arto Klami. Prior Knowledge Elicitation: The Past, Present, and Future. Bayesian Analysis, pp. 1
33, 2023.
12

---

Published as a conference paper at ICLR 2026
Petrus Mikkola, Luigi Acerbi, and Arto Klami. Preferential normalizing flows. In Advances in Neural Information Processing Systems, volume 37, pp. 55667-55702, 2024.
Kim A Nicoli, Christopher J Anders, Tobias Hartung, Karl Jansen, Pan Kessel, and Shinichi Nakajima. Detecting and mitigating mode-collapse for flow-based sampling of lattice field theories. Physical Review D, 108(11):1 14501, 2023.
Anthony O' O'Hagan. Expert knowledge elicitation: Subjective but scientific. The American Statistician, 73:69-81, 2019.
Long Ouyang, Jeffrey Wu, Xu Jiang, Diogo Almeida, Carroll Wainwright, Pamela Mishkin, Chong Zhang. Sandhini Agarwal, Katarina Slama, Alex Ray, et al. Training language models to follow instructions with human feedback. In Advances in Neural Information Processing Systems, volume 35, pp. 27730-27744 2022.
Art B. Owen. Monte Carlo theory, methods and examples. 2013.
R Kelley Pace and Ronald Barry. Sparse spatial autoregressions. Statistics & Probability Letters, 33 (3):291-297, 1997.
James Requeima, John Bronskill, Dami Choi, Richard Turner, and David Duvenaud. Llm processes: Numerical predictive distributions conditioned on natural language. In Advances in Neural Information Processing Systems, volume 37, pp. 109609-109671, 2024.
Danilo Rezende and Shakir Mohamed. Variational inference with normalizing flows. International Conference Machine learning, pp. 1530-1538, 2015.
Robin Rombach, Andreas Blattmann, Dominik Lorenz, Patrick Esser, and Björn Ommer. High- resolution image synthesis with latent diffusion models. In Proceedings of the IEEE/CVF Confer- ence on Computer Vision and Pattern Recognition, pp. 10684-10695, June 2022.
Murray Rosenblatt. Remarks multivariate transformation. The annals of mathematical statistics, 23(3):470472) 1952.
Anuj K Shah and Daniel M Oppenheimer. Heuristics made easy: an effort-reduction framework. Psychological Bulletin, 134(2):207, 2008
Yang Song and Stefano Ermon. Generative modeling by estimating gradients of the data distribution. In Advances in Neural Information Processing Systems, 2019.
Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. In International Conference on Learning Representations, 2021.
Vincent Stimper, Bernhard Schölkopf, and José Miguel Hernández- Resampling base distri- butions of normalizing flows. International Conference on Artificial Intelligence and Statistics, pp.
2022.
Richard S Sutton, Andrew G Barto, et al. Reinforcement learning: An introduction. 1998.
I. L. Thurstone. A law of comparative judgment. Psychological Review, 34(4):273-286, 1927.
Hugo Touvron, Louis Martin, Kevin Stone, Peter Albert, Amjad Almahairi, Yasmine Babaei, Nikolay Bashlykov, Soumya Batra, Prajjwal Bhargava, Shruti Bhosale, et al. Llama 2: Open foundation and fine-tuned chat models. arXiv preprini arXiv:2307.09288 2023.
Kenneth E Train. Discrete choice methods with simulation. Cambridge university press, 2009.
Pascal Vincent. A connection between matching and denoising autoencoders. Neural computa- tion, 23(7):1661-1674, 2011.
Bram Wallace, Meihua Dang, Rafael Rafailov, Linqi Zhou, Aaron Lou, Senthil Purushwalkam, Stefano Ermon, Caiming Xiong, Shafiq Joty, and Nikhil Naik. Diffusion model alignment using direct preference optimization. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pp. 8228-828888888 2238, 2024.
13

---

Published as a conference paper at ICLR 2026
Max Welling and Yee W Teh. Bayesian learning via stochastic gradient Langevin dynamics. In Proceedings of the 28th International Conference on Machine learning, pp. 681-688, 2011.
Kaiwen Zheng, Cheng Lu, Jianfei Chen, and Jun Zhu. Improved techniques for maximum likelihood estimation for diffusion odes. In International Conference Machine learning, pp. 42363-42389, 2023.
14

---

Published as a conference paper at ICLR 2026
APPENDIX
This appendix provides additional technical material complementing the main paper. Appendix A presents extended theoretical results for the exponential noise RUM. Appendix B contains the proofs. C Appendix provides method details and implementation specifications. Appendix D shows the RUM invariance under space reparameterization a uniform distribution. Appendix reportts experimental details, the runtime breakdown, and additional ablations on RUM model misspecification. Appendix F includes plots of the learned belief densities and a description of the LLM experiment. Appendix G describes the use of L.MMs in the preparation of this paper.
LIST OF NOTATION
x
vector
X
A matrix
Element of vector tt observation
X,
X
vector-valued random variable
W
RUM noise, a scalar random variable or a stochastic process
((x)
A stochastic process at index
A binary striet preference relation
x x'
x is chosen over x', where x is winner point and x' is loser point
p(x)
probability density function evaluated x
Belief density
A
Sampling density
Pw
Marginal winner density (MWD)
Px-xx
Joint density of winner-loser pairs encoded in the order of the arguments
u
Utility function
Tx(x)
The gradient of f with respect to x, and evaluated at Noise-conditioned
s0
network with parameters 0 Noise-conditioned score network evaluated at winner-loser pair (x, x') and noise scale
S x',
s
RUM noise level, positive scalar
r
Tempering field, a scalar field T X' - (0, oo) or a positive scalar (in case of constant field)
Txx)
Tempering field evaluated at x, a positive scalar
r*
Optimal tempering constant, positive scalar Weighting function, scalar field w X' - |0, oc)
(X)
Stochastic weight, scalar random variable l ---mm (Euclidean norm)
-1
15

---

Published as a conference paper at ICLR 2026
A TEMPERING FIELD UNDER THE EXPONENTIAL NOISE RUM
In this appendix, state the existence and provide closed-fome expression for the tempering field when the expert's choice model follows the exponential RUM with W ~ Exp(s).
Theorem A.1. Assume W ~ Exp(s). tempering field rxx) exists between the belief density p and the MWD pw, and it is given by the formula
2vol + 2x*(x) SU, "x''xxx
St,, pe (x')dxx :)
7(x)
(A.1)
(x)x) So p (xx)xx + +
where the sublevel set Lx = fx' E &t p(x) 2 p(x')3 and the superlevel set Ux = Lxx
Proof Sketch. For a fixed x E X, after lengthy manipulations, we get T(x) such that vx log pxx) r(x) Po Pw (x) 0. To handle the change in integration domains, we apply the generalized Leibniz integral rule (Flanders 1973). See Appendix B for the full proof.
L
6 for Bradley Terry and = for exponential RUM, both resulting in unit variance for case of comparison. The tempering fields extremely similar, but tempering in high-density regions appears slightly more pronounced in the Bradley- -Teryy model, due to a lighter-tailed conditional choice distribution.
Fig. A.1 illustrates the tempering fields from Theorem A.1 and Theorem 3.1 using
Figure A.1: Illustration of the tempering fields under two different RUMs when p Twomoons2D
(Stmmper et al. 2022). The tempering field r(x) of the exponential RUM (left, Theorem A.I and the Bradley-Tere Terry model (right, Theorem 3.1
B PROOFS
Theorem A.1 Assume W ~ Exp(s). A tempering field rxx) exists between the belief density p and the MWD Pan and is given by the formula,
(x)) Su. (xx'x' (x) St, (x' -1) 1)
2vol(Lx) + 2p 2p*(xx) (x) JUz "=x('="
7(x)
where the sublevel set Lx fx' X p(x) p(x')} and the superlevel set Ux (x' E X p(x) < p(x)))
Proof. We want to show that for each x E X, there exists a scalar r(x) > 0 such that log p(x)
T(x) V log Pw (x) 0. Our constructive proof determines this scalar through brute-force calculation. Fix a point x E X, and denote a constant rxx) > O.
Under the exponential RUM, the MWD P (x) equals
(x) 2x(x) / Fiuplace(0.1/)) (logp(x) logp(x')) X(x')dx'.
(B.1)
16

---

Published as a conference paper at ICLR 2026
For uniform X, this implies
/ xx
Vx log Px) Vx Vxlog/
Fiaplace(0.1/s)) (logp(x) log dx'.
(B.2)
Let 4x, Ux C X' be the regions 4x {x' E Ex p(x) 2 p(x)] and Ux - fx' E E 12(x) < <((x)}
Straightforward manipulations give us,
T > log Pro (x) - Vx log
(1ss
1 700(x)0 dxx)'' + f (1 (1-2) 2D 2 "xxx)=(xx²) dx')
1
= log
VV lxg
Ux
1xx
vol(Lx) + (x) Sü p- "xx')xxx - 20= (x) St., p² (x))xxx
= =vVx log
(xx
Because
Vx f(x) = Vx Txfxx)/f(x) the above vanishes, if and only if,
1 Vx (x)oo + 50 -0- (x) Ju
/ p "x(('' - - (x) /. pl(x' (x')dxx) O.
We will apply to LHS the generalized Leibniz integral rule (Flanders) 1973) for each fixed dimension j E [1....d] by interpreting ihe, as differentiation with respect to the time. To justify the generalized
Leibniz integral rule, note that the boundaries 0T 0Ux are defined by an implicit function
g : X' - X' whose graph is the set ((x, g())) = {(xx x') E X22 f(x, x') 0}, where the function x(x,x)== == p(x) p(x') is continuously differentiable by Assumption 3 Moreover, since fx x') = - Vp(x') * 0 almost everywhere, the implicit function theorem implics that the level
set ((x, x') x(xxx') = 0} locally defines x' as differentiable function of almost everywhere.
Therefore, (x) is continuously differentiable almost everywhere.
For smooth function f X R4 consider,
/ f(x')xx',
0x, A
where Ax = Lx Ax Ux. Interpret the scalar Xy time, and apply the generalized Leibniz integral rule, 2 0x, Ax
f(x))xx'= / f(x)dxx++ + f f(x'((n v)ds, A. ox,
Joa.
where n is the unit normal vector pointing outwards from the boundary 0Ag, v is the Eulerian velocity of the boundary aAx when x, is interpreted time, and dS is the surface element. The first term in RHS vanishes. For the second term, consider the level set [x' E X' p(x) - p(x)) = The normal vector is orthogonal to this level set, which equals to gradient with respect x = Vx (p(x) pxx'))/[[V (Vx))//VVV(((XXXX)X>>/>>>
p(x'))l = / V the normal vector of Lx while the normal vector of Ux is with the opposite sign.
This corresponds to
0x]
8x'y ), the velocity of the boundary with respect to xj. The total
Let us consider v =
(xx)'
... dx,
derivative of the boundary should not change,
Dx, (x(x) p(x))))))))) =0
d
Ox - x(x())
(p(x) - p(x'))
+ Σ (p(x)
0x,
Ox
0x,
i=1
0
Ox!
p(x)
plx
0.
0x,
Xx!
0x
Taking the dot product of the constraint with the normal vector gives,
0x, p(x)
n
Vx'((')I
17

---

Published as a conference paper at ICLR 2026
Since p(x') p(x) on 0Lx dUxx
0
8
(xx))xx==
-pp((x)
Öxx p(x)
1. || dS(x')
0x,
Lx
0
①
p
(xx (x) p(x)
JOL |[xx(x')] d5(x')
0x,
Ux
{ dx' = θ
0x, p(xx
V dS(x').
Oxy Lxx
JOLn
Together these imply that,
p 0x, lx
dx' + (x)
(x) exz Usx
/ p "x(('' =00" - a
I
(x())))))= = 0.
نxz lxx
We are left with the following terms that vanish,
p - (x) (-rvol(Lx)) +
1 p" (x) { p "x"(x'
2 =p=0(xx / 1'(x')dxx)
1
2 -
Ux
In order to this hold, the scalar in the parenthesis must vanish,
2 -vo((Lx) + -0((x) Ux
/ p "=((') - -xx) / lxx
p'(x')ctx'
= sxx(x) Use
p "=x('' + ss spp (xx) p'(x')
1.x
Rearranging the terms give us.
2vol (Lx) + 2p" (x) Stz p "xx(x" dxx T p'(x) "(x') + + +x=((x) Se, p px()) :)
Theorem 3.1
Assume W Gumbel(O, s). A tempering field rxx) exists between the belief density and the MWD Pw, and it is given by the formula,
Jx Itr, (x,x)) dxx
rx(x,x) (xx) = Sxx (1+r,(xx)) dx!
(B.3)
where (x,x') == på (x')p (x) is 1/ 1/sttempeeed density ratio.
Proof. Under the Bradley-Terry model, W ~ Gumbel(O, s), and the MWD (x) now equals to P= (x) = 2X(x) /
(oogiiste(o,s)) (log - ((og((x))
(B.4)
Following same lines of reasoning as in the constructive proof of Theorem A.1 we fix x E X and aim find constant 7(x) = > 0 solving the tempering field condition. Because p is uniform, the necessary and sufficient condition for the existence of r is that it solves the equation,
2 (x') 1 +(1 + 11+==((x))) (x))p (x))
(x')
/
p x (x)
dx' =0
This is equivalent to,
fx I+r.Ix.xT dx
==8
(B.5)
Sx T. xx(xx')) dxx
(1+7,
[p(x') 1/s pix)
where for clarity we denote xxx density ratio between density values at compared points x' and
pà (x')p (x)
which is 1/ s/s-tempered
18

---

Published as a conference paper at ICLR 2026
Proposition 3.2
Assume that there exists tempering field 7(x) between and q. A tempering defined by,
parameter T
Exap (x)X)=(X)).
(B.6)
where stochastic weight 2 O is given by
Vlogg(X))
w(X)
(I log (1)))
Ex-p
minimizes the Fisher divergence between p and q,
arg min F(p,
(B.7)
T>0
Proof. By the Leibniz formula,
0
V logg (x) - Vx log p(x)] p(x)dx
0
Vx logo (x) V log [x)]] p(x)dx
Ūr
/ 2
T V (x))
(V (x), Vx logp(x))) p(x)dx.
Since the Fisher score is quadratic in r, the critical point corresponds the global minimum, By the assumption, (Vx log q(x), V log p(x)) T(x) || V log (xx) Hence,
E Ex<p ([ V| log ((x)) Ex-p (IVlog((X)11))
Lemma B.4. Let T(x) be a tempering field. For any T >> holds that
F(p.q") / |- 1 V [x)] p(x)xx
Proof. Since 7(x) is tempering field,
V log (x) - Vx logg
|l(T -(x)) V lxg + (7(x)V>logg(X Vx lolgg((x))
[|(T - (x)) (x)]].
Proposition B.5. Let 7* be the optimal tempering and T(x) tempering field.
- E T(X) || log x(x))1²)
F(p.9" "")."")
E V log (x)))
19

---

Published as a conference paper at ICLR 2026
Proof. By Proposition 3.2 and Lemma B.4
Xx log 2" (x) - Vx log pl
p(x)dx
x
R
(IT T(X)12 1 2(
E [7(0x]]1000 -xx
2
E
E [Vlog((X'))1)
[IVlogg(X)Il
(T(X') [Vog
V log 0( (x10x((
E
E
(IVlog((X))1
E
((X) V (log((x)11))
E (X() ((x))//a log ((
+ E bog)
2
E
V log (x))()))1²²
I (lo (x))²
LE (T(X) V (o log
E >
E [IVlog
Proposition
3.3 Let r(x) be a tempering field. For any it holds that
F(p,4") = Ex-p (Ir (x-)((x)(= V (x))
(B.8)
Further, when T* > is the optimal tempering
(Exap (T(X) || log (og ((x))
F(p.g")) Ex-p (IVlog((x))
(B.9)
E Ex~p (I Vl q( (x))
Proof. This is combined result of Lemmal B.4 and Proposition B.5
Corollary B.7. Assume that the expert choice model follows the Bradley-Tery Terry model or the expo- nential RUM. The scores of the belief and the MWD are collinear. That is, there exists scalar-valued function rxx) such that,
V llog((x) T(x) V log (x).
Proof. The result follows directly from Definition 5|d Theorems 3.1 and A.1
C METHOD
C. TRAINING JOINT AND MARGINALS DISTRIBUTIONS
Recent works discuss in detail how diffusion models can be used to learn between joint and arbitrary conditional distributions, whereas modeling the marginals not always straightforward (Gloeckler al.) 2024). We adopt a simplified approach to estimate the marginal score function by leveraging corruption-basee marginalization strategy.
To model simultaneously both the joint distribution Pxxxx (x, x') and the marginal distribution Pwo (x), we introduce a binary conditioning variable joint E (true, false] into the score model. During training, randomly set joint false with 50% probability, and in this case, mask the input x, by replacing it with Gaussian noise (000, ofI), where o, is the current noise level. We then compute the denoising matching loss only over the winner dimensions (i.e., the first d
20

---

Published as a conference paper at ICLR 2026
(a) Tempered MWD, w/o joint training
(b) Tempered MWD, w/ joint training
Figure C.I: Replication Fig. using (a) only winner samples, with the score model trained only for the MWD (x). The quality of the final density estimate, i.e., the 'tempered" MWD, is clearly inferior compared to (b) training on the full joint Px-x (x, x') using both winners and losers.
Table C.I: Wasserstein distance between the winner marginal samples using sø (x, x', j,joit
0, temp = 0) and (x, x',
joint 1, temp 0). The last column reports the distance as
percentage relative to the fourth column of Table !
((x)
wasserstein
relative proportion
Onemoon2D
0.085 (#: 0.009)
23%
Ring2D
0.077 (# 0.012)
18%
Twomoons2D
0.072 (#: 0.010)
18%
MixturegaussiansSDD
0.071 (#: 0.002)
5%
Onegaussian
0.069 (#: 0.001)
6%
Stargaussian6D
0.156 (# 0.001)
12%
Mixturegaussians LOD
0.366 (#: 0.001)
28%
Onegaussian 161)
0.671 (+ 0.001)
13%
components of V log Px- x'). When joint true, train the model to predict the full joint score over both x and x'
At sampling time, to generate samples from the marginal distribution pw (x) using ALD, we similarly set joint = false se and replace x with Gaussian noise (O, If at each iteration of ALD. Fig, C.1 demonstrates the benefit of training the score model the full joint data while using the proposed marginalization method to model the marginal score.
Does the diffusion model learn to marginalize? To empirically study how well our approach for the diffusion marginalization is able to learn to marginalize, we analyze the difference between the marginal samples from the joint model, and the samples from the marginal model. Specifically, Table C.1 reports the Wasserstein distance between the samples from (i) the true winner marginal of the joint model: x where (x, x') Px-x (x,
(x,x') x') using the score model Sy (x, x', joint 1, temp 0), and (ii) the winner marginal of the joint model: x where
p=.(x) using the score model sy (x, x', J., joint = O, temp = 0). To quantify the similarity, we can contrast the Wasserstein distance between these two estimates to the one we have between the target and the score- score-based estimate (reported in the fourth column in Table [p. This relative difference is typically below 20%. indicating that the marginalization method performs well.
21

---

Published as a conference paper at ICLR 2026
C.2 DETAILS ON MODELING THE 'TEMPERED' MWD
To learn the MWD Pw (x), having only access to samples from the joint distribution of winners and losers Px-xx (x,x') x'), propose learning the full joint and the first marginal to capture the preference relationships in the data while enabling sampling from the tempered marginal via score-scaled ALD. To that end, we parametrize the score model So (x, x', jojont, temp) such that (x, x') + 8x (xxx x', join=11 temp=() models the score of the joint distribution of winners and losers, i.e. V log Px>xx/ (x, x'), and (x, temp) + Sx(x, 0, min joint=0, temp) models the score of the MWD and its 'tempered' version, i.c., log (x(x) and (xx V log Pw (x), respectively. This allows training the score network with both winners and losers, while still enabling sampling from the belief' density via the approximation (x 0, min, joint-O, temp=1) x r(x)V V log (x) V log p(x) by Corollary B.7 Fig. C. validates that the joint score model is superior to directly learning the marginal from only the winner samples.
We implement the method by training a score model through the denoising score matching equation on the concatenation of winners and losers, with random masking of losers during training (for more details, see Appendix C...n Technically speaking, to enable sampling via reverse diffusion and noise distribution (o) during training, defined
[
ALD, we a cosine noise schedule and a LogNormal( P Paeann P with mixture weight & = 0.5 assigned to the Dirac delta component. We stay as close as possible to the EDM-style diffusion model (Karras)
mixture of Dirac delta on
et al. 2024a). Specifically, we use the perturbation kernel Po x) = [xx I²I), which aligns with EDM and defines a forward diffusion process from mmin to max where pr (x) p(x) and Po (x) R (O, x I). After training Sx x', o, joint, 0) on perturbed winners and losers, we
sample from the belief density using the score-scaled ALD, We can either stop here, or optionally
train the tempered marginal score network SomnD (x, a, temp) whose weights can be initialized to that of so, through denoising score matching using the sampled data. Finally, we use the loss weighting
((()) o². Algorithms 2 summarize the method.
C.3 SCORE MODEL
We follow as closely as possible the HDM2 specifications used in the 2D toy experiment in (Karras
ltt al., 2024a). For both the joint score network and the tempered MWD score network, we use an MLP with one input layer and four hidden layers, SiLU activation functions (Hendrycks &
Gimpel 2016 are applied after each hidden layer, and implemented using the magnitude-preserving primitives from EDM2 (Karras et al. 2024b). In the joint score network, the input is (2d + 3)- dimensional vector (z, oo,.nt, temp),
and the output of each hidden layer has h features, where h E {32, 64, 96, 128} depending on the experiment dimensionality. In the MWD score network, the input is a (d +3) -d---miiisii vector (x. a, O. temp). The binary variables joint and temp are linearly embedded into h/ 8/4-dimensional space. Further, simple residual connection is applied to the embedded variables through all hidden layers. Otherwise, we use the same preconditioning for the score network as described in EDM (Karras et al. 2022).
C.4 BELIEF DENSITY RATIO MODEL
We parametrize the belief density ratio (x, x') RC p(x')/p(x) via parameterizing the unnormalized log density fx) Ri log p(x) + constant as an ML. with three hidden layers with SiLu activa- tions, and one output layer, such that log ro ro(xx xx)
fo(x') fx)).. The number of hidden units is tied to that of the score model (Section C.3) Regularization of the weights 0 is important for obtaining sensible results. To this end, we apply adaptive fz-regularization using the Adam opti- mizer (Kingma & Ba/ 2015) with weight decay. In contrast, standard & f--regularization, corresponding to Adam W (Loshchilov & Hutter 2019), yielded slightly inferior empirical performance. We set the weight decay to 10 except in the small data n = 100d experiments, where we use a higher value of 3 10-3
22

---

Published as a conference paper at ICLR 2026
C.5 TEMPERING HIELD ESTIMATE
We estimate the integral ratio in Eq. 6 by approximating both the numerator and the denominator using Monte Carlo estimates, using the MWD pio (x) as importance sampler,
Σ (x;) (x.))
-(x) R
(C.1)
Σ (xc) (x,
rolx.x:)
where the sums are over the importance samples Pur This Monte Carlo ratio estimator is biased, but consistent.
For the d-dimensional target, we use 2000d importance samples to estimate the integrals in the tempering field. The importance weights are computed using the probability-fie ODE of the MWD diffusion model (see Appendix for details). The estimated tempering field r(x) is clipped such that 1 S Txx) 2 Q- (0.99,, where Q- (0.99) denotes the 99% quantile of the estimated tempering field values (for visualization, use the 99.9% quantile). The lower bound follows directly from the theory (i.e., from formula (). while the upper bound is introduced for numerical stability to remove outliers.
C.6 SCORE-SCALED ALD
Score-sealed ALD uses the scaled score (x)) (x),, where V log is replaced by our estimated score. While 7(x) is not the ALD step size e 0, it is clear that rxx) influences the ALD update in manner similar to e. To ensure convergence of score-scaled ALD, at each ALD step
the step size (x) musx" where Cbase is the base step size be specified. The required number of iterations T should be chosen with respect to ébase- In our experiments, we use
o²
50, T = 40, and Cbuse = 0.15, except in the exeexerriments, where we use Caase 7.0 with I, = 15 and = 40. In the 2 we keep the original domain scale ([-3, 3) and do not rescale it to [-0.5, 0.5] unlike in the other experiments.
Regarding the injection of a deterministic ALD noise schedule during denoising score-matching training. find that the cosine schedule yields better empirical performance, while the noise schedule corresponding to the EDM time-step discretization is also a natural option.
C.7 DENSITY EVALUATION OF A DIFFUSION MODEL
Chen al. (0018) showed that for random variable whose probability density evolves over time, with dynamies dx f(x,, f) dt (where f is Lipschitz continuous in and continuous in f), the log-density at point Xo is given by
log (x = log pr (XT) + (- V f(xe, dt,
(C.2)
where V- denotes the divergence operator, which is equal to the trace of the Jacobian, In practice, the divergence often approximated using the Skilling Hutchhnson trace estimator (Grathwohl et al. 2019
Fex)===-oioo (x)))
where the expectation is typically estimated using a finite number of samples.
Now, applying this to EDM-type diffusion model, which characterized by the probability- flow ODE
dx o V log (x do,
(C.3)
we can compute the probability density at point x as
Fraex
logp(x) ~ log N (X=... 0, 0 m I)
V (x do,
(C.4)
where the score V ) is approximated by the score network So .), and the divergence term estimated using the Skilling- Hutchisson estimator.
23

---

Published as a conference paper at ICLR 2026
The ODE in Eq. C.3 is often stiff at small noise scales, requiring careful choice of numerical integration method for stable results. To integrate Eq. C.4 we solve the coupled system
dx - - (x,
log p(x) o sx(x,
do
do
where the latter ODE tracks the accumulation
log-density.
To solve the coupled system, we use the implicit Sllams-ashforth---aa implemented in the torchdiffeq package (Chen 2018). We compute the divergence exactly using Hutchinson estimator also worked, sometimes even yielding better results. Finally, to further stabilize the density estimates, we clamp the importance weights 1/p(x) to the interval defined by their 1st and 90th percentiles.
black-box ODE solver
automatic differentiation, although the
D RUM UNDER THE SPACE REPARAMETEERIATION
This appendix studies the exponential RUM and the Bradley-Terry model under space reparame- terization, and verifies that both RUMs are invariant under this transformation. This justifies our approach of transforming possible non-uniform X(x) into uniform distribution.
The invariance to the space reparamcterization holds for Fechnerian RUMs with the winner density,
(x) 21(x) / F (log log A(x))xxx
where F is a cumulative distribution function of the choice probability. Let a transformation T push the sampling density À into the uniform distribution on the hypercube, that Tz À Unif((0, 1Jd). It is enough to show that the winner density in the transformed space has the same form but with uniform sampling density and different belief density. To that end, denote y = (xx and note that
== x=XTT (y)) det (y)).
Hence, X(x) = Xy)- in the front of the integral cancels the volume change under the transfor- mation T. The integral in the new coordinate system y' can be computed by the change of variable Sx Gyd((T" 14) = Jo.yy [Gy Tdd for (x'x) == F (logp(T"> (x))""(()) log p(x')) and
denotes the Lebesgue measure on the hypercube [0, 1]d. The winner density in the transformed space reads as
Pw,trans 2
(log (logp(T-'(y)) - log loggppTT
dy'.
10.111
That is, the winner density in the transformed space follows the same RUM than in the original space, but with uniform sampling distribution and belief density Ptrans (y) pyT-l(y)). Note that the
normalization of Ptrane (y) is irrelevant, as the constant cancels out within the logarithmic utility difference. Consequently, the scale of the RUM noise (parameterized by F') remains invariant under the space reparameterization.
24

---

Published as a conference paper at ICLR 2026
E EXPERIMENTAL DETAILS
E. TARGET DISTRIBUTIONS
The log unnormalized densities of the target distributions used in the synthetic experiments are provided below.
| -
22
X1 +2 2
Onemoon2D
0.2
2
0.3
[||x||
1)2
(xx1
2)²
4x1 am
Twomoons2D
+ log 1 +
0.08
0.18
Ring2D log
Σ
32
12[ ))
where k =1
Stargaussian6D
1 log
NNx(1 E1)
2 N(x 11, E2,2)
a2
pa²
pa?
pa
g²
po
pa?
a² = 1, p 0.9, d = 6, 31d, EΙ =
po
Do
0²
po
po
pa
a²
pa?
pa
1) d-
pa²
poo
a²
pa
d po
2
po?
o²
d po
Σg
(-1)d-1 2a²
(-1)d-2po² po?
(-1)d-3pol 2
1 (x- mi) TE, (x (i)) m
Mixturegaussians, E {4,10)
log
exp
i=1
where
μ. v
3,
V1 1a,
V2 - 1a, [(-1)]] j=1 = V3.
Σ Q diag o ....,
Q o = 0² 0.1, Q [i | Rxd
Gaussian, DE (4,16)
-x n1) T Σ (x 11) 1 2
1)a
10
15
15
Σ
15
10
15
15
15
TO
E.2 OTHER EXPERIMENTAL DETAILS
Hyperparameters and optimization details. The score models are trained for varying number of iterations (d 8192, 2 < d < 10 : 12288, d 2 10 15360) with the Adam optimizer (Kingma Ba 2015) and batch size min(n, 4000) pairwise comparisons, where n the number of pairwise comparisons in the dataset. For the 2D experiments, we follow (Karras et al. 2024a) and use an adaptive learning rate, specifically a decay schedule of are/ max(iter, iter 1). with fnwer-function EMA profile with dre 0.01.
are/ 0.005 and iter,e/ / 1024 iterations. We The setup is somewhat sensitive to hyperparameters, and performance can vary depending on their tuning. We expect achieve better worse performance in the experiments depending how well
25

---

Published as a conference paper at ICLR 2026
Table E.1: Robustness to RUM noise family and noise level. The method score-r(x) assumes the Bradley-Terry model with noise level 0.7797, while the data mata-generating process is varied across three RUMs: Exponential RUM (W Exp(s)), Bradley-Terry (W (,mbel s)) and Thurstone- Mosteller (W (0, ))), each evaluated at noise levels E {0.1, 0.7797, 1.0, 5.0). Reported values are averages over 10 repetitions of the Wasserstein distance (4).
0.1
0.7797
1.0
5.0
Onemoon2D
W Exp(s)
0.693
(=0.078)
0.266
(+0.122)
0.252
(=0.122)
0.242
(=0.117)
W Gu s)
0.358
(=31.38)
0.366
(+5)145)
0.372 (=0114
0.599
+..132)
W Normal (0, )
0.294 (19.
0.318 (10.124)
0.321
(10.122)
0.549
(10.119)
Twomoons2D
W Exp(s)
0.577 (=0.069)
0.348
(25))
0.371
(=0.12
0.383
( 1339
W (0,
0.481
(=0.090)
0.440
(=0009
0,446 (=0.094)
0.500
00000)
W ~ Normal (0, s²)
0.375 (10.144)
0.345
(10.158)
0.346
(100160)
0.447
(..00)
Ring2D
W ~ Exp(s)
0.557 (=00071)
0.423
(0006)
0.427
(=006)
0,453
(=0 0.0.8
W G0, s)
0.462 (20)128)
0.387 (0.07))
0.368
(0.074)
0.449
(0.04)
W Normal )²)
0.495 (=0.10
0.430 (0.00)
0.420 (0.08)
0.491
(0.07)
the hyperparameters are tuned. The chosen hyperparameters are likely suboptimal, and we expect performance gains, especially in higher-dimensional experiments, if the hyperparameters are well tuned.
Environment. All experiments are conducted on server equipped with nodes containing dual Intel Xcon Cascade Lake processors (20 cores each, 2. IGH)). While exact training times and memory usage were not recorded, the datasets and score network architectures used are relatively lightweight.
Experiment replications. Every experiment was replicated with 10 different seeds, ranging from to 10.
Baseline. We used the official implementation of (Mikkola et al. 2024) and the provided config files to match the hyperparameter configuration used in their experiments to the closest experiment in paper. For example, for 2D experiments, we use the config file that was used in their Onemoon21) experiment.
E.3 ROBUSTNESS TO RUM NOISE FAMILY AND NOISE LEVEL
To study robustness of the method for misspecification in the data-generation process, we rerun Onemoon2D, Twomoons2D, and Ring2D by varying the true data-generation RUM noise family W and the noise level 8. In all cases, our method assumes the Bradley-Terry model with fixed noise level s = 6672. Table E.1 reports the results.
The results suggest that lower true noise levels generally lead to higher-quality estimates (note that for W ~ Exp(s), smaller corresponds to a higher noise level). However, using a true noise level that roughly matches the assumed model often yields better or comparable performance than choosing a very small noise level (note that the model is correctly specified only when W ~ Gumbel(O, s) This observation aligns with the theoretical identifiability of the problem discussed in Section 2.3 Overall, the results appear relatively robust to misspecification of both the noise family and the noise level.
E.4 RUNTIME BREAKDOWN
The computation times are profiled for the experiments in Section 5 Table E.2
reports the total
runtime allocated to three main components: DSM training of the score network (Section 4. :
estimation of the tempering field, including training the density ratio model re (Section 4.2 ; and
score-scaled ALD sampling (Section 4.3) with the number of samples varying from 25k to 40k. The final column shows the average time required to generate single sample. In practice, this
26

---

Published as a conference paper at ICLR 2026
also depends the sampling batch size, which we have not optimized and which is constrained by available memory.
Table F.2: Breakdown of computation times (in seconds). The table reports mean and standard deviation over 10 replicates.
sampling
DSM training
T(x) estimation
total
per sample
Onemoon2D
80 (= 16)
49 (# 12)
690 (# 45)
0,0276
Ring2D
82 (# 15)
47 (=: 10)
692 (= 41)
0.0277
Twomoons2D
( 14)
48 (+: 10)
701 42)
0.0280
Mixturegaaussians4D
262 (# 51)
4560 (+: 873)
6958 743)
0.2783
Onegaussian44D
283 (#: 46)
2660 27:)
7187 404)
0,2872
Stargaussian6D
440 (d: 54)
3364 496)
12053 (29)
0.4018
Mixturegaussians 1OD
653 (4:47)
6728 (4 1030)
27314 (4 853)
0.6829
Onegaussian16D
797 (#: 24)
4977 (# 223)
47477 (# 727)
1,1870
27

---

Published as a conference paper at ICLR 2026
F PLOTS
F.1 PLOTS OF LEARNED BELIEF DENSITIES
(a) Score-based
(b) Flow
Figure F.1: Onemoon2D experiment. The target distribution is shown heatmap, and samples from the learned model are overlaid in blue. (a) Samples from the score-base model. For this particular seed, the estimated tempering field is not too far from the true field, resulting in good fit. (b) Samples from the flow model.
(a) Score-based
(b) Flow
Figure F2: Twomoons2D experiment. The target distribution shown as a heatmap, and samples from the learned model are overlaid in blue. (a) Samples from the score-base model. For this particular seed, the estimated tempering field is not too far from the true field, resulting in good fit. (b) Samples from the flow model.
F.2 LLM EXPERIMENT
Details of the data generation for the LLM experiment, including the prompts used, are provided in Mikkola al. (2024. Appendix C.2). We reuse their data and scripts to convert the 5-wise rankings into the pairwise comparisons assumed in our setup:
Fig. F.10 completes the partial plot in main text for the LLM experiment. The elicited 2D and 1D marginals have the same support as the true data distribution marginals, and their shapes are also similar, with the distinction that score-based methods tend to generate Gaussian-like marginals. The only exception is the variable AveOccup, whose marginal appears to have an unreasonably long tail.
Fig. F compares our score-based diffusion method and the flow method of learning the LLM prior from pairwise comparisons in the LLM experiment, highlighting that the score-based method results in smoother estimates than the flow method. Table F. summarizes the densities quantitative manner by reporting the means for all variables.
28

---

Published as a conference paper at ICLR 2026
(a)
(b) Flow
Figure F.3: Ring2D experiment. The target distribution is shown as heatmap, and samples from the learned model are overlaid in blue. (a) Samples from the cooe-bbase model. (b) Samples from the flow model.
-5.0
25
59
-50
-2.
do
25
59
-50
-25
60
25
50
-35
25
5.
-35
an
35
50
-3.5
th
35
Figure F.4: Gaussian4D experiment. The target distribution is depicted by light blue contour points and its marginal by pink curve. The learned diffusion model depicted by greenish blue contour sample points and its marginal by a black curve.
29

---

Published as a conference paper at ICLR 2026
2-
2-
Figure F.5: Mixturegaussians4D experiment. The target distribution depicted by light blue contour points and its marginal by a pink curve. The learned diffusion model is depicted by greenish blue contour sample points and its marginal by a black curve.
30

---

Published as a conference paper at ICLR 2026
10
10
10
10
Figure F.6: Stargaussian6D experiment. The target distribution is depicted by light blue contour points and its marginal by a pink curve. The learned diffusion model is depicted by greenish blue contour sample points and its marginal by a black curve.
31

---

Published as a conference paper at ICLR 2026
10
Figure F.7: Stargaussian6D experiment when using the baseline flow-based method. The target distribution depicted by light blue contour points and its marginal by pink curve. The learned flow model is depicted by greenish blue contour sample points and its marginal by a black curve.
32

---

Published as a conference paper at ICLR 2026
Figure F8: Mixturegaussians 10D experiment. The target distribution is depicted by light blue contour points and its marginal by a pink curve. The learned diffusion model is depicted by greenish blue contour sample points and its marginal by a black curve.
33

---

Published as a conference paper at ICLR 2026
Figure F.9: Gaussian 16 16D experiment. The target distribution is depicted by light blue contour points and its marginal by a pink curve. The learned diffusion model depicted by greenish blue contour sample points and its marginal by black curve.
34

---

Published as a conference paper at ICLR 2026
Learned LLM from
Calroma detaset dientbution
Figure F.10: Full result plot for the LLM expert elicitation experiment, complementing the partial plot presented in Fig.
Learned
1-1
Figure F.11:
Comparison of the results when learning the LLM prior from pairwise comparisons using our score-based diffusion method (left) and the flow based method (right).
Table F.1: The means of the variables based (first row) the distribution of the California housing dataset, (second row) the sample from the score-based diffusion model fitted to the LLM's feedback, and (third row) the sample from the flow model.
MedIne HouseAge AveRooms AveBedrms Population AveOccup
Lat
Long
True data
3.87
28.64
5,43
1.1
1425.48
3,07
35.63 -119,57
Score-based
5.89
27.56
6,66
1.53
2997.96
4.70
36.69
-119.37
Flow
5.83
28.48
6.68
1.49
2948.17
3.36
36.73 -119.30
35

---

Published as a conference paper at ICLR 2026
G THE USE OF LARGE LANGUAGE MODELS (LLMS)
The data for Experiment "LLMM as a proxy for the expert" Section was obtained by prompting an LLM (Claude Haiku by Anthropic, March 2024). Further, the first version of the Rosenblatt transformation implementation was developed using an LLM. This version was later improved, and the final version was verified correctly transform points to the hypercube. The inverse transformation also worked in the tested cases. Finally, an LLM was used to check for writing and content errors in both the text and the code.
36


Reviews:
[{"review_id": "EtxWHpKcjl", "reviewer": "Reviewer_NgZe", "summary": "In this paper, the authors propose a density estimation method for recovering the target density from pairwise comparison data. The key contribution is showing that the score of the target density is colinear with the score of the winner’s density, allowing it to be estimated via score matching. The proposed method is evaluated on synthetic data generated from a Random Utility Model (RUM) as well as data derived from large language model (LLM) experts, demonstrating promising performance compared with Mikkola et al. (2024).", "questions": "What is the definition of the expert's belief density? why is it useful to generate sample from it?\n\n4.1, I feel the main estimation objective (i.e., the denoising score matching) should be displayed here, showing how the data distribution is constructed from winners and losers and how exactly they are used in the score matching objective.  The current explanations from line 324 - 328 and B1 is not clear enough. For example, what does \"train sθ(x,x′,σ,joint) using score-matching\" mean?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 2, "strengths": "The problem of estimating the density function from expert knowledge is novel in generative model context and underexplored. \n\nThe score factorization in (6) is a natural extension of Mikkola et al. (2024)'s proposal where tau(x) here is a constant. The connection with score matching (Fisher divergence) and generative model is very natural as well. \n\nThe performance reported in Section 5 is a clear improvement over the normalising flow model used in Mikkola et al. (2024). \n\nThe use of LLM to generate expert belief is an efficient way of producing cheap expert comparisons and is a motivating problem setting for the proposed method. \n\nThis paper is reasonably well-written (although I think the problem should be better motivated, see below).", "weaknesses": "The main weakness is the confusion of the setting considered in this paper. In what application, it is useful to estimate the experts beliefe density? For example, if I have experts opinions {x > x'} where x is a hand-written image and experts are asked copmare images based on how much they look like a digit one. In what application, would I want to estimate p(x)? I guess the user would probably be interested in knowing whether the experts consider a specific image x being 1 or not but I cannot see why it is useful to sample from p(x)? \n\nI could not find much motivation on that in the introduction section, and the experiment section is largely based on toy datasets. \n    - In the housing data example, it is unclear what is the \"belief about the distribution of features in the California housing dataset\". Does it mean, the expert belief on whether the data point is in the dataset or not? \n\nThat said, I can see that mathematically, this is an interesting object, as the method seems to estimate the utility function u(x) in RUM model (line 169), which is often unobserved. \n\nIn the experiments, the comparisons are made against Mikkola et al. (2024) using a flow model. However, as a benchmark, the authors should also include a diffusion model with \\tau(x) treated as a constant, in order to compare the performance gained by incorporating the estimated \\tau(x).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "oDfVgPVii9", "reviewer": "Reviewer_qpz9", "summary": "The paper studies density estimation from pairwise comparisons, linking the unobserved target density to a tempered winner density. By learning the winner’s score through score matching and applying a de-tempering step, the target density can be recovered. The authors prove a collinearity relation between the belief and winner scores, derive a tempering field, and propose an estimator under the Bradley–Terry model. A diffusion-based method using annealed Langevin dynamics demonstrates effective learning of complex belief densities from limited comparisons.", "questions": "Please see Weakness.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "1. Estimating the density function from pairwise comparison is novel.\n2. The motivation and structure of the work are clear.\n3. The theoretical results can show the benefits of the proposed method.", "weaknesses": "1. The notation is a little messy. Is there a notation list table to explain the meaning of these items?\n2. Maybe it is better to provide an analysis for the efficiency of the proposed method. There should be some experiment results or some discussions on it.\n3. I recommend to provide a more clear structure of the algorithm in Section 4, which should be the core contribution of your work.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "tl2g2qWX5E", "reviewer": "Reviewer_aZ9i", "summary": "This paper studies density estimation from pairwise comparisons: For example, an expert is shown two i.i.d. draws $\\mathbf{x}, \\mathbf{x}'$ from a sampling distribution $\\lambda(\\mathbf{x})$, and the object they prefer is recorded. It is assumed that $\\mathbf{x}$ is preferred over $\\mathbf{x}'$ with probability $\\mathbb{P}(\\mathbf{x} \\succ \\mathbf{x}')$, which is modeled as a random utility model (RUM) with deterministic utility $u(\\mathbf{x}) = \\log p(\\mathbf{x})$. The goal of this work is to infer $p$ from a set of $(\\mathbf{x}, \\mathbf{x}', \\mathbf{x} \\succ \\mathbf{x}')$ triples. Firstly, for two special cases of RUM (Bradley-Terry & exponential noise RUM), the authors (1) provide theoretical results that relate the marginal winner density (MWD) $p_w(\\mathbf{x}) \\propto \\int \\mathbb{P}(\\mathbf{x} \\succ \\mathbf{x}') \\lambda(\\mathbf{x}) \\lambda(\\mathbf{x}') \\, \\mathrm{d}\\mathbf{x}'$ to $p(\\mathbf{x})$ via a position-dependent *tempering field* $\\tau(\\mathbf{x})$ by showing that their scores are collinear, i.e., $\\nabla \\log p(\\mathbf{x}) = \\tau(\\mathbf{x}) \\nabla \\log p_w(\\mathbf{x})$, and (2) provide explicit formulas for $\\tau(\\mathbf{x})$. \nSecondly, the paper builds on these insights by proposing a practical algorithm: The authors train (1) a diffusion model to estimate $\\nabla \\log p_w$ (MWD model), and (2) a second neural network to estimate the density ratio $r_\\theta(\\mathbf{x}, \\mathbf{x}') \\approx p(\\mathbf{x}')/p(\\mathbf{x})$, which is then used to estimate $\\tau(\\mathbf{x})$ (with the MWD model as a proposal distribution).\nLeveraging the relation to the score of $p$ allows the authors to then draw approximate samples from $p$ using score-scaled annealed Langevin dynamics. \nEmpirical findings show that this approach outperforms a previous flow-based method qualitatively and quantitatively.", "questions": "+ L373 claims that a \"linear MLP score network\" is used. What exactly is meant with \"linear\" here? The used MLP is surely a non-linear function of inputs/parameters.\n+ Have the authors considered training a model to output $\\tau(\\mathbf{x})$ directly, e.g. by regressing many (expensive) Monte Carlo estimates? This would amortize the computational complexity during sampling.\n+ In low dimensions (e.g., $d = 2$), numerical integration algorithms (e.g. quadrature) have much better convergence rates than Monte Carlo. Have the authors tried such methods to estimate $\\tau(\\mathbf{x})$ (instead of importance-sampled Monte Carlo)?\n+ When training the diffusion model on the full joint (Figure C.1 (b)), is the distribution of the model when `joint = False` actually close to the *true* marginal of the joint model? An experiment would be interesting where (1) you sample $(\\mathbf{x}, \\mathbf{x}')$ from the model with `joint = True`and discard $\\mathbf{x}'$ (which is an *exact* sample from the winner marginal), and (2) you sample $\\mathbf{x}$ with `joint = False`, and compare the distributions (e.g. Wasserstein).\n+ Why is it that the diffusion model trained on the full joint (both winners and losers) learns \"better\" winner marginals than the model only trained on winner? I find this surprising.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "+ **Motivation.** The work is well-motivated and tackles an important issue in preference learning. The studied RUM models are relevant in many disciplines.\n+ **Theoretical Results.** The theoretical results are sound, well-justified, novel, and important for this line of research. Moreover, they give rise to novel practical algorithms for estimating $p$ using score-based approaches.\n+ **Presentation.** The main results are presented in a clear way. It was easy to follow the author's main arguments.", "weaknesses": "+ **Estimating $\\tau(\\mathbf{x})$** \n\t+ Given a parametric model $r_{\\theta}(\\mathbf{x}, \\mathbf{x}') \\approx p(\\mathbf{x}')/p(\\mathbf{x})$, estimating $\\tau(\\mathbf{x})$ still involves approximating two (possibly high-dimensional) integrals (Eq. 7) for each $\\mathbf{x}$, which is computationally very expensive.  \n\t+ While estimating these integrals with importance-sampled Monte-Carlo (MC) yields unbiased estimates for each of the integrals, the resulting MC estimate for $\\tau(\\mathbf{x})$ is *not* unbiased, since it involves a fraction of the two integrals. The paper would benefit from a more thorough discussion of this estimation procedure. Moreover, claiming that \"the integrals are computed using importance sampling [...]\" (L343-L344) is misleading (\"computed\" should be replaced).\n\t\n+ **Computational cost for sampling from the estimated $p$**\n\t+ In score-scaled annealed Langevin dynamics (ALD), using the scaled score $\\tau(\\mathbf{x}) \\nabla \\log p_w(\\mathbf{x})$ requires estimating $\\tau(\\mathbf{x})$ in each step of ALD. In Appendix C.6, it is stated that for a $d$-dimensional target, $2000d$ importance samples are used to estimate the integrals in Eq. 7. Hence, a single ALD step is at least $2000d$ times more expensive than regular ALD sampling (where the score is obtained by a single forward pass of the network)---while not counting the computational complexity of calling $r_\\theta$ and estimating $\\log p_w(\\mathbf{x})$ with the probability-flow ODE.\n\t+ It is clear that this does not scale well to large $d$, or large networks (neither to large diffusion model, nor large $r_\\theta$).\n\t\n+ **Score collinearity when $\\sigma > 0$**\n\t+ The identity $\\nabla \\log p(\\mathbf{x}) = \\tau(\\mathbf{x}) \\nabla \\log p_w(\\mathbf{x})$ only holds in the *noise-free* case, i.e., $\\sigma = 0$. When replacing $p$, $p_w$ with their Gaussian-smoothed version $p*\\mathcal N(0,\\sigma^2I)$ and $p_w*\\mathcal N(0,\\sigma^2I)$ for some $\\sigma > 0$, this is *not* true anymore. I think this should be stressed more explicitly in the main text.\n\t+ L353-L354 claims that the \"tempering field relation (Eq. 6) is guaranteed to hold in the small-noise limit [...]\". This is misleading, as Eq. 6 only holds when $\\sigma = 0$ exactly, and is only approximate when $\\sigma \\approx 0$.\n\t+ L354: \"Possible mismatch at higher noise scales is not an issue [...]\". This is only true in the ideal case of $L \\to \\infty$. In practical settings, where $L$ is moderate, the scores at higher noise levels will have a large effect on ALD.\n\t\n+ **Hyperparameters** \n\t+ Many of the empirical results seem sensitive to the particular choice of both training and sampling hyperparameters. In Appendix C.7, it seems that $\\epsilon_{\\text{base}}$  is chosen very differently for 2D experiments ($\\epsilon_{\\text{base}} = 7.0$) than in the other experiments ($\\epsilon_{\\text{base}} = 0.15$). This could hint at the fact that the magnitude of $\\tau(\\mathbf{x})$ may be miscalibrated.\n\t\n+ L109: \"The network is trained to predict the score of this kernel, which is typically tractable.\"\n\t+ This is not true: When minimizing the objective in Eq. (1), the network will learn $s_{\\theta}(\\mathbf{x}, \\sigma) \\approx \\nabla_{\\mathbf{x}} \\log p_\\sigma(\\mathbf{x})$, which is the score of $p_\\sigma$, and not of the perturbation kernel $p_\\sigma(\\tilde{\\mathbf{x}} \\mid \\mathbf{x})$.\n\t\n+ Very minor suggestions regarding notation: \n\t+ L163: I suppose $C$ should be a set, so this should read $C = \\\\{\\mathbf{x}, \\mathbf{x}' \\\\}$ instead. \n\t+ L374-375 reads \"2$d$ and 1$d$\" marginals, but $d$ refers to the *particular* dimension of the target (L367), so it should rather read e.g. \"2D and 1D\" (where \"D\" is used as a shorthand for \"dimensional\").\n\t\n+ Typos:\n\t+ L1022, L1052: \"probability ODE\", should probably read \"probability-flow ODE\" \n\t+ L1049: The equation for the Skilling-Hutchinson trace estimator should read $=$, not $\\approx$. The expectation is the exact trace, estimating the expectation with Monte Carlo makes it approximate.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "gOTdnaO92e", "reviewer": "Reviewer_zzGP", "summary": "This paper provides a method to sample from a distribution given only independent paired samples from an arbitrary distribution and a ranking of pairs of points namely $(x,x',I[U(x)>U(x')])$  where $U$ is some *random* utility function defined as $U(x)=\\log p(x)+W$ with $W$ a known noise distribution. It achieves this in two stages: firstly, it learns the score marginal distribution of the sample with highest utility ($p_w(x)$) via denoised score matching, second it learns the transformation $\\tau(x)=\\frac{\\nabla\\log p(x)}{\\nabla\\log p_{w}(x)}$ in order to obtain an estimate for the data distributions score. The paper derives a closed form for $\\tau$ as a function of the density ratio when $W$ follows a Gumbel or Exponential distribution as a function of the belief density ratio, $p(x)/p(x')$ for any $x,x'$ pair.", "questions": "* I am unclear of the relationship between 3.1 and the rest of the paper. I appreciate that when $\\tau$ is constant you then get the relationship you establish later on ($\\nabla \\log p=\\tau\\nabla\\log q$) however the later work does not work under the assumption that $\\tau$ is constant.\n\n* Can you use the score at the end to produce an explicit density estimate or does it exclusively give a score function and thereby a means of sampling from a distribution.\n\n* Is there any way to use the score estimate to give a new proposal distribution say $\\lambda'$ which is closer to $p$?\n\n* Were other choices for estimation of $s_\\theta$ explored such as normalising flows?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 2, "presentation": 4, "contribution": 3, "strengths": "*  Well written with clear introductions to each of methods used.\n\n* Tackles an interesting and challenging problem where only pairwise ranks can be used to learn a multivariate distribution.\n\n* Good mix of synthetic illustrative experiments and real world experiments outside of the simplifying assumptions of the paper.\n\n* Provides theoretical justification explanation of  existing method for density estimation with ranked comparisons giving theoretical lower bounds on accuracy of prior work.", "weaknesses": "* Estimation of the transformation $\\tau$ is non-trivial and relies upon estimating the belief density ratio $r(x,x')=p(x)/p(x')$. Once the belief density ratio has been estimated you already have access to the unnormalised density (for arbitrary fixed $x'$ $p(x)=C\\cdot r(x,x'))$ for all $x$) and the score of the distribution ($s(x)=\\nabla\\log p(x)=\\lim_{h\\rightarrow 0} \\log r(x,x+h)/h$). As such it is not clear that the requirement to learn the MWD and then the tempering field $\\tau$ from $r$ is strictly necessary.\n\n* Estimation of $r$ and by consequence $\\tau$ is dependent upon the distribution of $W$ being known and from one of two distributional families. This should be made clear in the paper earlier stages of the paper. (It might also be good to highlight that the problem is unidentifiable if $x\\succ x'$ is taken as the deterministic density rankings as this is fundamental to the approach and without this the equation on line 79 doesn't seem to make sense as $x\\succ x'$ is deterministic.)\n\n* It is unclear how the quality of the method depends upon the accuracy of the assumption on $W$, both in terms of the parameter $s$ and the assumption of the distributional family (namely Gumbel). It would be good to have experiments exploring synthetic data in a more reasonable miss-specified setting.\n\n* The score matching using both winning and paired samples seems slightly under-explored, I appreciate the nice illustration of the improvement from including the paired samples but is there any further justification as opposed to only utilising samples from the distribution you are learning. It seems to me that for this to hold we would reasonably need the neural network to learn that the masking variable represents the marginalisation of the distribution in order for paired samples to meaningfully help learn the MWD.\n\n* (This is more a comment) Algorithms B.1 and B.2 should be in the main body if possible and clarification that the likelihood on line 342 comes from the assumption that $W\\~ \\text{Gumbel}(0,s)$.\n\n* Small mistake I believe $\\mathbb{P}(x\\succ x')=F_{W(x')-W(x)}(u(x)-u(x'))$ not $F_{W(x)-W(x')}(u(x)-u(x'))$ (although I could be mistaken.) It's also maybe worth highlighting that this probability is what you mean by the choice distribution conditional on $\\mathcal{C}$.\n\n\nHighlight early on that $x\\succ x'$ is random with some fixed noise distribution known up to a parametric family, as this is necessary for the identifiability of the problem.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)