REVIEW 2 major objections 6 minor 28 references
MCBench: A Benchmark Suite for Monte Carlo Sampling Algorithms
T0 review · 2 major / 6 minor · reviewed 2026-08-10 · deepseek-v4-flash
Pith's one-line read MCBench is a Julia benchmark suite that claims to make Monte Carlo sampling quality quantitatively comparable by pitting user samples against IID reference samples from known targets using sliced Wasserstein distance and maximum mean…
desk verdict Useful modular benchmark suite, but the default ESS estimator miscalibrates the comparison bands for MCMC output, which is the main use case. read the letter →
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
The reading
What carries the argument
The load-bearing mechanism is the batch-comparison protocol. Each sample set is divided into $m$ batches of effective sample size $n$; the default effective sample size estimator is the ratio of the squared sum of the weights to the sum of squared weights, which is intended to make autocorrelated sampler output comparable to $n$ IID samples. For each batch, the suite computes basic metrics (marginal mean, variance, chi-square) and, for two-sample comparisons, the sliced Wasserstein distance and the maximum mean discrepancy. The sliced Wasserstein distance projects both distributions onto random directions and averages the one-dimensional Wasserstein distances; the maximum mean discrepancy compares mean embeddings in a Gaussian reproducing kernel Hilbert space, computed exactly or via random Fourier features. The distributions of metric values over the IID batches define the expected scatter and the normalization bands used to judge the user-sampler batches.
What would settle it
Re-run the 3D Gaussian walkthrough twice: once with the default weight-based effective sample size estimator and once with an autocorrelation-aware estimator, for example one based on integrated autocorrelation time. If the Metropolis-Hastings metric markers shift substantially relative to the 1-, 2-, and 3-sigma IID bands, the suite's normalized comparison depends on the choice of ESS estimator rather than on sampling quality alone. A second check: run a sampler that provably misses one mode of the correlated mixture target and verify that the SWD and MMD markers fall outside the 3-sigma bands; if they do not, the metrics are not sensitive enough to detect the failure the paper claims they expose.
Extended reading notes
Core claim
The central claim is that MCBench provides a general, domain-neutral workflow for evaluating Monte Carlo sampling algorithms. The user picks a target from a list of IID-sampleable distributions (including correlated and multimodal Gaussians, Cauchy, and the eight-schools hierarchical model), generates samples with the sampler of their choice, and imports them into the suite. MCBench splits both the user samples and its own IID reference samples into $m$ batches of effective sample size $n$, evaluates each metric on every batch, and builds the reference distribution of metric values from the IID batches. The user-sampler batches are then normalized relative to that IID distribution, so the final plots show whether the sampler's metric mean sits inside the green, yellow, or red bands. The paper demonstrates the workflow with the Metropolis-Hastings sampler from BAT.jl: on an uncorrelated 3D Gaussian the metrics match the IID expectation, while on a correlated mixture model the metric distributions clearly deviate, showing the suite can expose a sampler's failure to capture multimodal structure.
Load-bearing premise
The load-bearing premise is that the default effective sample size estimator correctly splits a sampler's output into batches that behave like independent samples; the paper itself notes in Section 6.1 that this estimator overestimates the effective sample size for the Metropolis-Hastings chain, so the comparison bands could be narrower than they should be and the verdict could flatter a correlated sampler.
Editorial extensions
If this is right
- Users can benchmark external samplers by writing samples to a CSV or HDF5 file; no integration with MCBench's sampling code is required.
- The normalized bands give a thresholded visual test: a sampler whose metric mean falls outside the 3-sigma region fails to reproduce the target according to that metric.
- Because every target is IID-sampleable, any detected discrepancy is attributable to the sampler rather than to reference-sample noise.
- The modular TestCase-Sampler-Metric design means new target functions and new distance metrics can be added without changing the comparison protocol.
- The eight-schools target provides a realistic hierarchical model with three prior levels, letting users test samplers on a nontrivial posterior where an IID accept-reject reference is still available.
Reading between the lines
- The IID-vs-IID metric distribution functions as an empirical null distribution, so one could assign an empirical p-value to a user sampler's metric value; the paper displays bands but does not compute p-values.
- Replacing the default effective sample size estimator with an autocorrelation-aware one would likely widen the user-sampler bands for MCMC output, changing how close to the IID bands a sampler appears; this is a direct, testable consequence of the limitation admitted in Section 6.1.
- The same machinery could rank approximate inference methods beyond MCMC, such as variational approximations or normalizing flows, as long as they emit samples; the paper targets MC sampling but the comparison protocol does not depend on how the user sample was produced.
- MMD with a Gaussian kernel is a kernel two-sample test, so MCBench could be extended from distance reporting to calibrated power analysis, answering 'which metric detects a given sampler failure at which sample size?' rather than only 'how far is this sampler from IID?'
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. The paper presents MCBench, a Julia package for benchmarking Monte Carlo sampling algorithms. Users select one of several IID-sampleable target distributions, draw samples with an external sampler, and the package compares those samples against IID samples from the same target using a set of metrics: basic statistics (marginal mean, marginal variance, chi-square), the sliced Wasserstein distance (SWD), and the maximum mean discrepancy (MMD), optionally with random Fourier features. The comparison is visualized by normalizing the IID-batch metric distributions to mean zero and coloring the 1-sigma, 2-sigma, and 3-sigma regions, then overlaying the user sampler's metric mean and spread. Two walkthrough examples with the BAT.jl Metropolis-Hastings sampler are given: a 3D standard normal and a mixture of correlated normals. The appendix contains algorithms and theoretical background for SWD and MMD.
Significance. If correctly calibrated, MCBench would fill a practical need for a domain-neutral, modular, and easy-to-use benchmark suite for Monte Carlo samplers. The paper provides an open-source Julia package, clear workflow diagrams, example implementations, and reproducible walkthroughs, which are valuable strengths. The theoretical appendices for SWD and MMD are standard and appear correct in their main claims. However, the default effective-sample-size estimator is miscalibrated for autocorrelated samplers, and the paper explicitly admits this in Section 6.1. Since the central purpose of the benchmark is to provide quantitative, reliable comparisons of sampling quality, this issue is load-bearing: the comparison bands are built from IID batches of nominal size n, while batches from MCMC samplers have a smaller true effective size, so the method can produce false negatives even for samplers that are exact at stationarity. The paper should be revised to correct this calibration problem before the quantitative comparison claims are supported.
major comments (2)
- [§2, §6.1] The default effective-sample-size estimator defined in Section 2 as the ratio of the squared sum of the sample weights divided by the sum of their squares is the Kish formula; for unweighted MCMC output it returns the raw sample count, independent of autocorrelation. Section 6.1 states that the observed overestimation of the metric variance is due to an overestimation of the effective sample size for the Metropolis-Hastings chain. Consequently, the user-sample batches are not of effective size n, while the IID reference bands are computed from batches of n independent samples. Even a sampler that is exact in stationarity will therefore have per-batch metric distributions that are wider than the colored bands, producing systematic false negatives for every autocorrelated sampler evaluated with the default settings. This directly undermines the central claim of enabling reliable, quantitative measures of sampling quality. The fix is to use an autocorrelation-aware ESS estimator (e.g., batch means or a spectral-variance estimator) as the default, to require the user to supply a valid ESS, or to construct both the IID and user-sample batches using the same effective-size definition; a validation study demonstrating the calibration on samplers with known ESS would also be needed.
- [§2 (two-sample metrics)] The workflow for the advanced two-sample metrics (SWD and MMD) is described in Section 2 as comparing IID batches directly with user-generated samples, without specifying how the effective sample size is handled. SWD and MMD estimates depend on the number of samples and on autocorrelation: a correlated sampler with the same raw count as an IID batch carries less independent information, so its metric values will have higher variance. If the user samples are not thinned or weighted to the same effective size as the IID reference batches, the comparison is miscalibrated in the same way as for the basic metrics. The manuscript should specify the effective-sample-size treatment for two-sample metrics, or explicitly state that only raw sample counts are matched and calibrate the bands accordingly.
minor comments (6)
- [§4.2.1 heading] The heading 'Sliced W assertein Distance' contains a typo; it should be 'Sliced Wasserstein Distance'.
- [§6.2] The code defines a 10-dimensional mixture (MvNormal(r*ones(10), ones(10,10)*0.9 + I(10)*0.1) and names it 'normal_3d_multimodal_10std'), but the section title and text say 'in 3D'. Please clarify the dimension used in the walkthrough.
- [§6.1 code] The parameters n=100 and n_samples=10^5 are not clearly defined; the text interprets them as 100 batches of 10^5 samples each, but the parameter names suggest n might be the batch size and n_samples the total number of samples. Please document the parameter semantics in Section 5.
- [Appendix B.2] In the RFF Monte Carlo approximation, the formula 'k(x,y) ≈ 1/D Σ 2 cos(ω^T + b_i)' is missing the feature arguments; it should read cos(ω_i^T x + b_i) cos(ω_i^T y + b_i).
- [Appendix B.2] The spectral density of the Gaussian kernel is stated inconsistently: Algorithm 3 uses ω ∼ N(0, σ^{-2}I), while the text first writes p(ω)=N(0,2σ^{-2}I) and then gives a density with exponent −∥ω∥^2 σ^2/2, which matches neither. Please unify the Fourier-domain normalization.
- [Table 1] The Eight Schools entry shows the likelihood but not the prior distributions from the text; including the full model in the table or caption would make the table self-contained.
Circularity Check
No circular derivation: MCBench compares user samples to externally generated IID baselines; the main flagged issue is the unvalidated default ESS estimator for autocorrelated samplers, which is a correctness caveat rather than a circular step.
full rationale
The paper makes no fitted prediction and derives no scientific constant. The benchmark's core is an external comparison: user-supplied samples are partitioned into batches and compared, via standard metrics (mean, variance, SWD, MMD), to batches of IID samples generated from the same target function. None of these metrics is defined in terms of the benchmark output, and the target distributions are standard or externally defined (posteriorDB eight schools). The default ESS estimator in Section 2 is stated as a default formula (squared sum of weights over sum of squares), not fitted to the metric results, so the user-batch construction is not circular even though Section 6.1 admits it overestimates ESS for Metropolis-Hastings and thereby widens the per-batch metric variance. The BAT.jl citation in Section 6 is an example sampler, not a load-bearing justification; the benchmark accepts FileBasedSampler data and works with Stan and PyMC as well. Because none of the paper's claims reduces to its own outputs or to an author-imposed uniqueness or ansatz, there is no significant circularity. The ESS issue belongs in a correctness or robustness review, not in a circularity finding.
Assumptions & free parameters
free parameters (3)
- MMD Gaussian kernel bandwidth sigma =
median of pairwise distances, data-dependent
- SWD number of projections L =
not fixed in the paper; appendix suggests about 50 projections are sufficient
- RFF feature dimension D =
not fixed in the paper; appendix uses D from 10 to 1000
assumptions (4)
- standard math Fourier slice theorem and invertibility of the Fourier transform on Rd
- standard math Bochner's theorem for continuous shift-invariant positive definite kernels
- domain assumption IID reference samples from the bundled target distributions are exact draws from those distributions
- domain assumption The default effective sample size estimate is adequate for batching correlated samples
Cite this review
Pith. "Pith review of MCBench: A Benchmark Suite for Monte Carlo Sampling Algorithms." pith.science (2026). https://pith.science/paper/5GEIQ4FR
@misc{pith2026250103138,
author = {Pith},
title = {Pith review of: MCBench: A Benchmark Suite for Monte Carlo Sampling Algorithms},
year = {2026},
howpublished = {\url{https://pith.science/paper/5GEIQ4FR}},
note = {Machine review of arXiv:2501.03138}
}
read the original abstract
In this paper, we present MCBench, a benchmark suite designed to assess the quality of Monte Carlo (MC) samples. The benchmark suite enables quantitative comparisons of samples by applying different metrics, including basic statistical metrics as well as more complex measures, in particular the sliced Wasserstein distance and the maximum mean discrepancy. We apply these metrics to point clouds of both independent and identically distributed (IID) samples and correlated samples generated by MC techniques, such as Markov Chain Monte Carlo or Nested Sampling. Through repeated comparisons, we evaluate test statistics of the metrics, allowing to evaluate the quality of the MC sampling algorithms. Our benchmark suite offers a variety of target functions with different complexities and dimensionalities, providing a versatile platform for testing the capabilities of sampling algorithms. Implemented as a Julia package, MCBench enables users to easily select test cases and metrics from the provided collections, which can be extended as needed. Users can run external sampling algorithms of their choice on these test functions and input the resulting samples to obtain detailed metrics that quantify the quality of their samples compared to the IID samples generated by our package. This approach yields clear, quantitative measures of sampling quality and allows for informed decisions about the effectiveness of different sampling methods. By offering such a standardized method for evaluating MC sampling quality, our benchmark suite provides researchers and practitioners from many scientific fields, such as the natural sciences, engineering, or the social sciences with a valuable tool for developing, validating and refining sampling algorithms.
Reference graph
Works this paper leans on
-
[1]
MIT Press, Cambridge, MA (2006)
Williams, C.K.I., Rasmussen, C.E.: Gaussian Processes for Machine Learning. MIT Press, Cambridge, MA (2006)
work page 2006
-
[2]
Meegen, C., Waschinky, N., Ickstadt, K., Ricken, T., Herbrandt, S., Henning, C.: Uncertainty quantification for fluid saturated porous media - Bayesian analysis, variational sensitivity analysis, surrogate modeling and reliability analysis. In: Kaliske, M., Graf, W. (eds.) Polymorphic Uncertainty Modelling for the Numeri- cal Design of Structures. Lecture...
work page 2025
-
[3]
The Journal of Chemical Physics 21(6), 1087–1092 (1953) https://doi.org/10.1063/1.1699114
Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H., Teller, E.: Equation of State Calculations by Fast Computing Machines. The Journal of Chemical Physics 21(6), 1087–1092 (1953) https://doi.org/10.1063/1.1699114
-
[4]
Journal of Statistical Software 76(1) (2017)
Carpenter, B., Gelman, A., Hoffman, M.D., Lee, D., Goodrich, B., Betancourt, M., Brubaker, M., Guo, J., Li, P., Riddell, A.: Stan: A probabilistic programming language. Journal of Statistical Software 76(1) (2017)
work page 2017
-
[5]
SN Computer Science 2(3), 210 (2021) https://doi.org/10
Schulz, O., Beaujean, F., Caldwell, A., Grunwald, C., Hafych, V., Kr¨ oninger, K., Cagnina, S.L., R¨ ohrig, L., Shtembari, L.: BAT.jl: A Julia-based tool for 16 Bayesian inference. SN Computer Science 2(3), 210 (2021) https://doi.org/10. 1007/s42979-021-00626-4
work page 2021
-
[6]
Grossi, S., Letizia, M., Torre, R.: Refereeing the Referees: Evaluating Two-Sample Tests for Validating Generators in Precision Sciences (2024) arXiv:2409.16336
arXiv 2024
-
[7]
arXiv preprint arXiv:2407.04967 (2024)
Magnusson, M., Torgander, J., B¨ urkner, P.-C., Zhang, L., Carpenter, B., Vehtari, A.: posteriordb: Testing, benchmarking and developing Bayesian inference algo- rithms. arXiv preprint arXiv:2407.04967 (2024)
arXiv 2024
-
[8]
Monge, G.: M´ emoire sur la th´ eorie des d´ eblais et des remblais. Mem. Math. Phys. Acad. Royale Sci., 666–704 (1781)
Show all 28 references
-
[9]
In: Dokl
Kantorovich, L.V.: On the translocation of masses. In: Dokl. Akad. Nauk. USSR (NS), vol. 37, pp. 199–201 (1942)
1942
-
[10]
In: International Conference on Machine Learning, pp
Arjovsky, M., Chintala, S., Bottou, L.: Wasserstein generative adversarial net- works. In: International Conference on Machine Learning, pp. 214–223 (2017)
2017
-
[11]
Advances in Data Analysis and Classification (2024) https://doi.org/10.1007/s11634-024-00599-1
Ding, Z., Omlor, S., Ickstadt, K., Munteanu, A.: Scalable Bayesian p-generalized probit and logistic regression. Advances in Data Analysis and Classification (2024) https://doi.org/10.1007/s11634-024-00599-1
2024 doi
-
[12]
Annual review of statistics and its application 6(1), 405–431 (2019)
Panaretos, V.M., Zemel, Y.: Statistical aspects of Wasserstein distances. Annual review of statistics and its application 6(1), 405–431 (2019)
2019
-
[13]
Advances in Neural Information Processing Systems 35, 36985–36998 (2022)
Nguyen, K., Ho, N.: Amortized projection optimization for sliced Wasserstein generative models. Advances in Neural Information Processing Systems 35, 36985–36998 (2022)
2022
-
[14]
In: Scale Space and Variational Methods in Com- puter Vision: Third International Conference, SSVM 2011, Ein-Gedi, Israel, May 29–June 2, 2011, Revised Selected Papers 3, pp
Rabin, J., Peyr´ e, G., Delon, J., Bernot, M.: Wasserstein barycenter and its application to texture mixing. In: Scale Space and Variational Methods in Com- puter Vision: Third International Conference, SSVM 2011, Ein-Gedi, Israel, May 29–June 2, 2011, Revised Selected Papers ...
2012
-
[15]
The Journal of Machine Learning Research 13(1), 723–773 (2012)
Gretton, A., Borgwardt, K.M., Rasch, M.J., Sch¨ olkopf, B., Smola, A.: A ker- nel two-sample test. The Journal of Machine Learning Research 13(1), 723–773 (2012)
2012
-
[16]
The Journal of Machine Learning Research 11, 1517–1561 (2010)
Sriperumbudur, B.K., Gretton, A., Fukumizu, K., Sch¨ olkopf, B., Lanckriet, G.R.: Hilbert space embeddings and metrics on probability measures. The Journal of Machine Learning Research 11, 1517–1561 (2010)
2010
-
[17]
Advances in neural information processing systems 25 (2012) 17
Gretton, A., Sejdinovic, D., Strathmann, H., Balakrishnan, S., Pontil, M., Fuku- mizu, K., Sriperumbudur, B.K.: Optimal kernel choice for large-scale two-sample tests. Advances in neural information processing systems 25 (2012) 17
2012
-
[18]
Bochner, S., Functions, M., Integrals, S., Analysis, H., Tenenbaum, M., Pol- lard, H.: Lectures on Fourier Integrals. (AM-42). Princeton University Press, ??? (1959). http://www.jstor.org/stable/j.ctt1b9s09r
1959
-
[19]
Journal of Statistical Software 98(16), 1–30 (2021) https://doi.org/10.18637/jss.v098.i16
Besan¸ con, M., Papamarkou, T., Anthoff, D., Arslan, A., Byrne, S., Lin, D., Pearson, J.: Distributions.jl: Definition and modeling of probability distributions in the juliastats ecosystem. Journal of Statistical Software 98(16), 1–30 (2021) https://doi.org/10.18637/jss.v098.i16
2021 doi
-
[20]
https://doi.org/10.5281/ zenodo.2647458
Lin, D., White, J.M., Byrne, S., Bates, D., Noack, A., Pearson, J., Arslan, A., Squire, K., Anthoff, D., Papamarkou, T., Besan¸ con, M., Drugowitsch, J., Schauer, M., contributors: JuliaStats/Distributions.jl: a Julia package for prob- ability distributions and associated func...
2019 doi
-
[21]
Journal of Open Source Software 8(87), 5236 (2023) https://doi.org/10.21105/joss.05236
Roualdes, E.A., Ward, B., Carpenter, B., Seyboldt, A., Axen, S.D.: Bridgestan: Efficient in-memory access to the methods of a stan model. Journal of Open Source Software 8(87), 5236 (2023) https://doi.org/10.21105/joss.05236
2023 doi
-
[22]
Schuhmacher, D., B¨ ahre, B., Bonneel, N., Gottschlich, C., Hartmann, V., Heine- mann, F., Schmitzer, B., Schrieber, J.: transport: Computation of Optimal Transport Plans and Wasserstein Distances. (2024). R package version 0.15-4. https://cran.r-project.org/package=transport
2024
-
[23]
Prentice Hall, Upper Saddle River, NJ (2000)
Munkres, J.R.: Topology, 2nd edn. Prentice Hall, Upper Saddle River, NJ (2000)
2000
-
[24]
PhD thesis, Universit´ e Paris Sud-Paris XI; Scuola normale superiore (Pise, Italie) (2013)
Bonnotte, N.: Unidimensional and evolution methods for optimal transportation. PhD thesis, Universit´ e Paris Sud-Paris XI; Scuola normale superiore (Pise, Italie) (2013)
2013
-
[25]
In: Uncertainty in Artificial Intelligence, pp
Xie, Y., Wang, X., Wang, R., Zha, H.: A fast proximal point method for computing exact Wasserstein distance. In: Uncertainty in Artificial Intelligence, pp. 433–453 (2020)
2020
-
[26]
In: ICASSP 2024-2024 IEEE Inter- national Conference on Acoustics, Speech and Signal Processing (ICASSP), pp
Le, D., Nguyen, H., Nguyen, K., Nguyen, T., Ho, N.: Fast approximation of the generalized sliced-Wasserstein distance. In: ICASSP 2024-2024 IEEE Inter- national Conference on Acoustics, Speech and Signal Processing (ICASSP), pp. 6920–6924 (2024)
2024
-
[27]
Advances in neural information processing systems 19 (2006)
Gretton, A., Borgwardt, K., Rasch, M., Sch¨ olkopf, B., Smola, A.: A kernel method for the two-sample-problem. Advances in neural information processing systems 19 (2006)
2006
-
[28]
Bodenham, D.A., Kawahara, Y.: eummd: efficiently computing the mmd two- sample test statistic for univariate data. Statistics and Computing 33(5), 110 (2023) 18 Appendix A Algorithms Algorithm 1 Calculating the sliced Wasserstein distance Require: Two sample sets: X = {xj}N j=...
2023
Reviewed August 10, 2026 · model on record in the stance chip above.
Discussion (0). Continue with ORCID to comment.