Approximating $f$-Divergences with Rank Statistics

AI summary

Overview

• Research area: Statistical machine learning and information theory; specifically divergence estimation between probability distributions, high-dimensional probability, and generative modeling.
• Technical level: Advanced. The paper is theory-forward, with proofs of monotonicity, lower-bound, convergence rates, finite-sample deviation bounds, and asymptotic normality, followed by empirical benchmarks and a generative transport demonstration.
• Scope: The paper proposes and analyzes an optimization-free, order-statistics-based estimator of f-divergences built from rank histograms, together with a sliced high-dimensional extension, empirical benchmarking against neural baselines, and use as a learning objective.

What This Paper Is About

Standard f-divergences such as Kullback-Leibler, total variation, Hellinger and chi-squared are defined through the density ratio dμ/dν, which is hard to estimate reliably in moderate-to-high dimensions. The authors take a different route: they compare distributions using only order information by discretizing the probability integral transform into a "rank histogram" with K bins and measuring how far that histogram is from uniform through a discrete f-divergence. The resulting estimator requires no density-ratio estimation and no iterative optimization, and is then lifted to R^d by averaging over random one-dimensional projections.

Key Contributions

1. A rank-histogram approximation of f-divergences. The authors define (D^{(K)}_{f,\nu}(\mu)), a resolution-K approximation that is an optimization-free estimator of the (sliced) f-divergence, generalizing prior work by de Frutos et al. (2024a, 2024b, 2026) to arbitrary entropy functions (f), including differentiable ones that interact well with automatic differentiation.
2. Theory for the one-dimensional construction. They prove (D^{(K)}{f,\nu}(\mu)) is nondecreasing in (K), bounded above by the true (D{f,\nu}(\mu)), convex (and weakly lower semicontinuous when (R_\nu) is continuous), consistent as (K \to \infty) under mild regularity of the quantile-domain density ratio, with quantitative rates; they also prove finite-sample deviation bounds and asymptotic normality.
3. Sliced rank-statistic f-divergences in R^d. They define the sliced version by averaging the univariate construction over random unit-direction projections, proving the chain (D^{(K)}{f,\nu}(\mu) \le SD{f,\nu}(\mu) \le D_{f,\nu}(\mu)), convexity, and convergence of the sliced approximation to the sliced divergence.
4. Empirical validation. They benchmark against classical and neural baselines on synthetic tasks, showing stable approximation in high dimensions and in small-sample regimes, and demonstrate use as a sample-based transport objective for implicit models on two-dimensional toy problems and on image datasets (CIFAR-10).

Main Findings

• Monotonicity and lower-bound property. Theorem 2.3 establishes (D^{(K)}{f,\nu}(\mu) \le D^{(K+1)}{f,\nu}(\mu) \le D_{f,\nu}(\mu)) exactly. The construction is convex, and weakly lower semicontinuous when (R_\nu) is continuous. A Markov-kernel interpretation via Bernstein polynomials ((b_{n,K})) gives an alternative route through the data processing inequality for f-divergences.
• Connection to prior work. With the total-variation generator (f_{TV} = |\cdot - 1|), Definition 2.2 recovers (up to a prefactor) the ISL discrepancy from de Frutos et al. (2024a): (D^{(K)}{f{TV},\nu}(\mu) = (K+1), d_K(\mu,\nu)), and the monotonicity result generalizes de Frutos et al. (2024a, Thm. 2) and de Frutos et al. (2026, Thm. 2.2).
• Convergence rates depend on smoothness of the quantile-domain density ratio. Theorem 2.5 gives: convergence to zero if (r \in C([0,1])); (O(K^{-\alpha/2})) if (r \in C^{0,\alpha}([0,1])); and the fast rate (O(K^{-1})) if either (r \in C^2([0,1])), or (r) is Lipschitz and (f \in C^2([0,\infty))).
• Which divergences get the fast rate. Example 2.2 notes (f_{TV}) is the only globally Lipschitz entropy function among the standard ones considered. The (O(K^{-\alpha/2})) rate applies to most divergences (including KL, Jensen–Shannon, squared Hellinger, and Jeffreys) when (0 \notin \mathrm{ran}(r)), since they are Lipschitz away from zero. The (O(K^{-1})) rate excludes KL and Hellinger and (|\cdot - 1|^\alpha) for (\alpha \in (1,2)), but holds for the chi-squared divergence with (f_{\chi^2}(t) = \frac{1}{2}(t-1)^2), other polynomial ("Tsallis") entropy functions, and the triangular discrimination generator (f(t) = \frac{(t-1)^2}{t+1}) (Lindsay, 1994).
• Explicit finite-sample guarantees. Theorem 2.6 gives, for fixed (K) and (f) that is (L_f)-Lipschitz on ([0, K+1]), an expected error bound (L_f (K+1)\sqrt{2\pi}\left(\frac{1}{\sqrt{N}} + \frac{1}{\sqrt{M}}\right)), and a concentration bound (L_f (K+1)\sqrt{2\log(2/\delta)\left(\frac{1}{N} + \frac{1}{M}\right)}) holding with probability at least (1-\delta).
• Sliced extension inherits the theory. Theorem 3.2 shows convexity of the sliced map and (D^{(K)}{f,\nu}(\mu) \le SD{f,\nu}(\mu) \le D_{f,\nu}(\mu)) for all (K \in \mathbb{N}), with (\lim_{K\to\infty} D^{(K)}{f,\nu}(\mu) = SD{f,\nu}(\mu)) when the projected quantile-domain density ratios are continuous for almost all directions.
• Asymptotic normality and an oracle choice of K. Theorem 3.3 gives (\sqrt{N}\left(D^{(K)}{f,\nu}(\hat{\mu}N) - D^{(K)}{f,\nu}(\mu)\right) \to \mathcal{N}(0, \tau_K^2)) in distribution when (f \in C^1([0,K+1])) with Lipschitz (f'). Theorem 3.4 shows that, under (C^2) regularity of (r) and bounded (f''), choosing (K_N = N^\beta) with (1/2 < \beta < 1) (equivalently (\sqrt{N}/K_N \to 0) and (K_N/N \to 0)) yields asymptotic normality with variance (\mathrm{Var}\mu\big(f'(r(R_\nu(X)))\big)).
• Competitive against a neural KL estimator. In the dimension benchmark replicating the synthetic setup of the neural KL estimator of Sreekumar et al. (2021) — Gaussian versus uniform, with (d \in {2, 5, 10}), sample sizes (10^4), (10^5), (10^6) and (5 \cdot 10^6), (K = 64), averaged over 10 independent runs with (\pm 1) standard deviation bands — the rank-statistic estimator (using an axis-corrected sum of per-coordinate terms) tracks the analytic KL reference closely and becomes increasingly stable as sample size grows, notably for (d = 5) and (d = 10). The neural estimator shows larger variability and more noticeable deviations from the reference, especially in higher dimensions.
• Ratio-to-reference behavior in 1D. Table 1 reports (D^{(K)}{f,\nu}(\mu)/D{f,\nu}(\mu)) for (K \in {32, 64, 128, 256, 512}) with (n_\mu = n_\nu = 10{,}000) and 10 seeds. Accuracy improves monotonically with (K) in every row. The hardest case is the heavy-tailed JS row (Laplace vs. Gaussian), with ratios 0.488, 0.651, 0.778, 0.869, and 0.933 across the five K values. Ratios slightly above one (up to 1.068 in the multimodal JS, (\Delta = 0.5) row at (K = 512)) are attributed to finite-sample Monte Carlo and numerical-reference error.
• Sliced ratios stay near one. With (K = 64), (L = 128), (n_\mu = n_\nu = 10{,}000), and (R = 10) runs, the ratio (d \cdot D^{(K)}{f,\nu}(\mu)/D{f,\nu}(\mu)) stays close to one with moderate variability across dimensions 2, 5, 10, 20, and 50, indicating that the (d)-scaling gives a reasonable normalization. Deviations grow more noticeable in higher dimensions, especially for JS and KL, while squared Hellinger remains comparatively stable.
• Trade-offs in choosing K and L. Small (K) gives a coarse but statistically stable rank histogram; large (K) gives a finer approximation at the cost of finite-sample variability, mirroring the classical binning trade-off in calibration-error estimation (Roelofs et al., 2022). For fixed (K), the projection error from using (L) directions behaves like (O(L^{-1/2})), as in general sliced probability divergences (Nadjahi et al., 2020); localized or highly anisotropic discrepancies may need substantially more projections, where structured or variance-reduced directions may be preferable.

Methodology in Plain English

The construction rests on the probability integral transform (Rosenblatt, 1952; Savage, 1952; Knothe, 1957): if (X \sim \nu) and (F_\nu) is the CDF, then (F_\nu(X)) is uniform on ([0,1]). Uniformity diagnostics built on this idea — rank histograms — are standard tools in forecast calibration and reliability assessment (Hamill, 2001; Gneiting et al., 2007).

The authors turn that diagnostic into a divergence. For a resolution parameter (K), the rank statistic (A^{(K)}{\mu|\nu}) is the number of (K) i.i.d. draws from (\nu) that fall at or below an independent draw from (\mu), taking values in ({0,\dots,K}). Its probability mass function (Q^{(K)}{\mu|\nu}) is a discrete "rank histogram": if (\mu = \nu) the histogram is uniform, and departures from uniformity signal mismatch. The rank-statistic f-divergence (D^{(K)}{f,\nu}(\mu)) is then the discrete f-divergence between this histogram and the uniform law (U_K), computed as (\frac{1}{K+1}\sum{n=0}^{K} f\big((K+1)Q^{(K)}_{\mu|\nu}(n)\big)). The whole thing reduces to sorting and counting, so no density-ratio model is needed and there is no inner optimization loop.

To prove the estimator behaves, the authors compare the histogram against a piecewise-constant version (r_K) of the quantile-domain density ratio (r), and show (r_K \to r) at a Bernstein-approximation rate, which transfers directly to the divergence gap. For high dimensions, they average the univariate divergence over random unit-direction projections on (\mathbb{S}^{d-1}) to get a sliced version, and they prove the same structural properties carry over.

For the generative experiments, they reinterpret the divergence as an energy. Particles are projected onto random directions; each projection's rank coordinates are computed with a smoothed empirical CDF that replaces indicator functions with logistic soft indicators; a Bernstein-smoothed histogram (\hat{Q}^{(K)}(\mathbf{v})(n) = \frac{1}{N}\sum_{i=1}^N b_{n,K}(v_i)) is formed; and the rank coordinates are refined by a proximal step that minimizes the divergence plus a trust-region penalty (\frac{1}{2\eta}|\mathbf{v} - \mathbf{v}_0|_2^2), approximated in practice by SGD or Langevin-type inner samplers such as ULA or MALA.

Baselines include neural variational f-divergence estimation, the Donsker–Varadhan representation used by MINE (Belghazi et al., 2018), and noise-contrastive estimation (Gutmann and Hyvärinen, 2010); the paper motivates the rank approach partly by the instabilities associated with adversarial and variational training (Arjovsky et al., 2017; Gulrajani et al., 2017) and the statistical error of density-ratio plug-in approaches (Moon and Hero, 2014; Rubenstein et al., 2019).

Why This Matters

The paper offers a divergence estimator that is bounded, depends only on order information, requires no density-ratio network, and is fully determined by the samples once (K) is fixed. That combination targets the main practical failure mode of f-divergence estimation: the need to fit flexible models of the ratio, which is statistically expensive in high dimensions and unstable under adversarial training. Because the construction is differentiable in suitable choices of (f), it also serves directly as a learning objective rather than only as a diagnostic.

Real-world applications touched by the paper:

• Generative modeling and implicit models: the divergence is used as a sample-based transport objective, demonstrated on two-dimensional toy problems and on image data (CIFAR-10; Krizhevsky, 2009).
• Forecast calibration and reliability assessment: rank histograms and probability-integral-transform uniformity checks are standard in this domain (Hamill, 2001; Gneiting et al., 2007), and this work turns that diagnostic into a general divergence family.
• Hypothesis testing and model comparison: f-divergences are described as arising throughout statistics, from hypothesis testing to model comparison (Ali and Silvey, 1966; Csiszár, 1967).
• Mutual-information and ratio-based estimation pipelines: the paper positions the method as an alternative to variational and neural estimators such as MINE and noise-contrastive estimation for settings where those are unreliable.

Industry relevance: the estimator avoids an iterative optimization loop at estimation time, scales to high dimensions through slicing, and remains usable in the small-sample regime where the neural baseline in the authors' benchmark had not yet stabilized — properties that matter for training objectives and monitoring metrics where compute and sample budgets are limited.

Future Directions

• Choosing projections better in high dimensions. The paper notes that a fixed (L) can lead to mild under- or over-estimation as (d) grows for JS and KL, and suggests structured or variance-reduced directions such as orthogonal or quasi-Monte Carlo projections, without evaluating them.
• Schedules for the resolution parameter. The authors propose progressive schedules that start from a coarse (K) and increase it as the coarse rank histogram stabilizes during training, but do not provide a general practical rule beyond the asymptotic oracle choice of Theorem 3.4.
• Solving the proximal rank-space step. Equation (9) has to be approximated in practice by deterministic updates such as SGD or by Langevin-type inner samplers such as ULA or MALA; the relative merits and convergence behavior of these inner solvers are left open.
• Extending beyond the current empirical scope. The benchmarking is on synthetic Gaussian/uniform, mean-shift, scale-change, Gaussian-mixture, and Laplace-vs-Gaussian families plus generative transport on two-dimensional toys and image datasets; whether the guarantees translate into advantages on other real distribution-matching tasks is not established here.

Target Audience

Researchers and practitioners in statistical machine learning, information theory, and generative modeling who work on divergence estimation, density-ratio-free methods, or sliced probability divergences. It is also relevant to statisticians interested in rank statistics, probability integral transforms, and calibration diagnostics, and to readers who want a theoretically grounded, optimization-free alternative to neural variational divergence estimators with explicit convergence rates and finite-sample bounds. A background in f-divergences and basic empirical process arguments will help; the construction itself is readable without it.

Authors’ abstract

We introduce a rank-statistic approximation of $f$-divergences that avoids explicit density-ratio estimation by working directly with the distribution of ranks. For a resolution parameter $K$, we map the mismatch between two univariate distributions $μ$ and $ν$ to a rank histogram on $\{ 0, \ldots, K\}$ and measure its deviation from uniformity via a discrete $f$-divergence, yielding a rank-statistic divergence estimator. We prove that the resulting estimator of the divergence is monotone in $K$, is always a lower bound of the true $f$-divergence, and we establish quantitative convergence rates for $K\to\infty$ under mild regularity of the quantile-domain density ratio. To handle high-dimensional data, we define the sliced rank-statistic $f$-divergence by averaging the univariate construction over random projections, and we provide convergence results for the sliced limit as well. We also derive finite-sample deviation bounds along with asymptotic normality results for the estimator. Finally, we empirically validate the approach by benchmarking against neural baselines and illustrating its use as a learning objective in generative modeling experiments.