Measuring Sample Quality with Stein's Method

1 Introduction

• Biased MCMC procedures can reduce Monte Carlo variance through faster sampling, but standard quality measures presume eventual convergence and therefore ignore asymptotic bias.This creates challenges for selecting samplers and their parameters.Source: §1 · p. 1
• The proposed quality measure compares sample and target expectations over a test-function class whose target expectations are known through Stein’s method.This avoids generic integration under the target distribution.Source: §1 · p. 1, §2 · p. 2
• The Stein discrepancy is computed by solving a linear program and is designed to compare asymptotically exact, asymptotically biased, and deterministic sample sequences.The paper also develops an efficient geometric-spanner implementation.Source: §1 · p. 1, §3 · p. 3
• The Stein operator depends on the target only through the derivative of its log density, so its normalizing constant need not be known.The operator is therefore computable when the target density is known only up to normalization.Source: §3.1 · p. 3
• For strongly log-concave densities with bounded continuous third and fourth derivatives, Stein discrepancy convergence to zero implies Wasserstein convergence.The stated class includes Bayesian logistic regression posteriors under Gaussian priors.Source: §3.2 · p. 3
• The classical Stein discrepancy converges to zero when sample measures converge in mean square and their target log-density gradients converge in mean.This provides a sufficient condition for convergence of the quality measure.Source: §3.3 · p. 4

Graph Stein Discrepancies

• The computable Stein discrepancy is formulated as an optimization program over functions g ∈G.Its objective is linear in g and depends only on function and gradient values at the discrete sample points.Source: §4.1 · p. 4
• The complete-graph Stein discrepancy is equivalent to the classical discrepancy up to a dimension- and norm-dependent constant κd.Specifically, S(Q, TP, G∥·∥) ≤ S(Q, TP, G∥·∥,Q,G1) ≤ κd S(Q, TP, G∥·∥).Source: §4.1 · pp. 4–5

Geometric Spanners

• Because the complete graph imposes O(n^2) constraints, sparse geometric spanners provide an O(n)-constraint alternative with comparable discrepancy values.For a t-spanner, S(Q, TP, G∥·∥,Q,G1) ≤ S(Q, TP, G∥·∥,Q,Gt) ≤ 2t^2 S(Q, TP, G∥·∥,Q,G1).Source: §4.1 · p. 5, §4.2 · p. 5
• A 2-spanner with O(κdn) edges can be computed in O(κdn log(n)) expected time, motivating its adoption as the standard quality measure.The constant κd depends only on dimension d and the selected norm.Source: §4.2 · p. 5

Implementation

• The recommended ℓ1 norm decomposes the multivariate optimization into d independent finite-dimensional linear programs that can be solved in parallel.For coordinate-boundary domains, additional compatibility constraints enforce smooth boundary conditions while retaining equivalence to the classical discrepancy.Source: §4.3 · p. 5, §4.4 · p. 5
• Algorithm 1 computes the multivariate measure using a sparse 2-spanner, parallel coordinate linear programs, and boundary constraints, while univariate cases use sorting instead.In one dimension, the complete-graph discrepancy can be computed by enforcing constraints only between consecutive sorted sample and boundary points.Source: §4.4 · pp. 5–6

Comparing Discrepancies

• For a Gaussian target, i.i.d. Gaussian samples have complete graph Stein discrepancy decay at n^-0.52, while scaled Student’s t samples remain bounded away from zero.The Student’s t distribution matches the Gaussian variance but has 10 degrees of freedom.Source: §5.1 · p. 6
• The classical Stein discrepancy can be computed in one dimension as a convex quadratically constrained quadratic program with O(n) variables and O(n) constraints.Source: §5.2 · p. 6

Selecting Sampler Hyperparameters

• ESS selects a larger SGLD step size than Stein discrepancy, with ESS maximized at ϵ = 5 × 10^-2 and Stein discrepancy minimized at ϵ = 5 × 10^-3.The experiment uses 50 independent sequences of length 1000 and selects the step size by median quality.Source: §5.3 · p. 7, Figure 3 · p. 8

Quantifying a Bias-Variance Trade-off

• Approximate RWMH has higher Stein quality under smaller computational budgets but is eventually overtaken by the asymptotically exact sequence.The comparison uses Bayesian logistic regression and 10^5 likelihood evaluations.Source: §5.4 · p. 7, Figure 4 · p. 8
• Additional normalized probability, mean, and second-moment error measures agree with the Stein discrepancy quantification in the approximate-sampling comparison.Source: §5.4 · p. 8

Assessing Convergence Rates

• For uniform samples, estimated Stein discrepancy rates are n^-0.49 for i.i.d. sequences, n^-1 for Sobol sequences, and n^-0.96 for kernel herding sequences.The i.i.d. and Sobol rates accord with expected O(1/√n) and O(1/n) bounds, while the herding rate exceeds its best known bound.Source: §5.5 · p. 8

Related diagnostics

• The proposed quality measure compares biased, exact, and deterministic sample sequences by exploiting an infinite class of known target functionals.Source: §6 · p. 8
• Finite-functional diagnostics can lose discriminating power; for a N(0, 1) target, the score statistic cannot distinguish samples with equal first and second moments.Source: §6 · p. 8

Deterministic sequences

• Figure 5 compares complete graph Stein discrepancy convergence for samples targeting P = Unif(0, 1).Source: Figure 5 · p. 9

Function distance

• The smooth function distance is defined using classes of real-valued functions with continuous derivatives.Source: B Proof of Theorem 2: Stein… · p. 9

Norms

• The general norm framework includes explicitly defined operator norms for vectors, matrices, and tensors.Source: B Proof of Theorem 2: Stein… · p. 9

Stein factors

• For strongly log-concave densities, explicit positive Stein factors support bounds relating the smooth function distance to the Stein discrepancy.Source: B Proof of Theorem 2: Stein… · p. 9

Stein equation

• The Langevin diffusion construction supplies a function solving the Stein equation for sufficiently smooth test functions with bounded derivatives.Source: B Proof of Theorem 2: Stein… · pp. 9–10

Wasserstein convergence

• The proof combines the Stein bound with a smooth-Wasserstein inequality to conclude Wasserstein convergence when the relevant Stein discrepancy converges to zero.Source: B Proof of Theorem 2: Stein… · p. 10

Expectation comparison

• For any admissible Stein function, the target expectation of its Stein operator is zero, enabling comparison between sample and target expectations.Source: C Proof of Proposition 3: Stein… · p. 10
• The expectation difference is decomposed into score-function and gradient terms involving differences between sample and target variables.Source: C Proof of Proposition 3: Stein… · p. 10

Norm bounds

• The second inequality follows from matrix norms, dual norms, Fenchel-Young, and Cauchy-Schwarz inequalities.Source: C Proof of Proposition 3: Stein… · p. 10

Conclusion

• Uniform validity over the admissible function class completes the proof.Source: C Proof of Proposition 3: Stein… · p. 10

• Scaling the test-function class by any positive constant does not change the Stein discrepancy.The objective is linear in g, so scaling the class preserves the normalized discrepancy.Source: D Proof of Proposition 4: Equivalence… · p. 11
• Function-class inclusion gives the first inequality relating the constrained and broader Stein discrepancies.The proof starts from G∥·∥ ⊆ G∥·∥,Q,G1.Source: E Proof of Proposition 5: Equivalence… · p. 11
• Whitney-Glaeser extension preserves function values and gradients on the sample support while introducing only the dimension-dependent factor κd.The extension theorem supplies a function in κdG∥·∥ with matching values and gradients at every sample point.Source: E Proof of Proposition 5: Equivalence… · p. 11
• Because the Stein objective depends only on sample-point values and gradients, the extension yields S(Q, TP, G∥·∥,Q,G1) ≤ κd S(Q, TP, G∥·∥).The preserved sample-point information makes the objective value unchanged before applying the scaling relation.Source: E Proof of Proposition 5: Equivalence… · p. 11
• The t-spanner property connects arbitrary support points through edges whose lengths control accumulated function and gradient variation.The proof applies the triangle inequality and tensor-norm definition along the spanner path.Source: F Proof of Proposition 6: Equivalence… · p. 11
• The spanner argument gives S(Q, TP, G∥·∥,Q,Gt) ≤ 2t^2 S(Q, TP, G∥·∥,Q,G1).Thus the discrepancy over the t-controlled class is bounded by a 2t^2 factor times the unit-scale discrepancy.Source: F Proof of Proposition 6: Equivalence… · p. 11

• Theorem 9 formulates the non-uniform classical Stein discrepancy on a sorted set containing sample and finite boundary points.The domain is X = (α, β), with −∞ ≤ α < β ≤ ∞ and sorted values x(1) < ··· < x(n′).Source: G Finite-dimensional Classical Stein… · pp. 11–12
• The program optimizes a Stein objective involving log p(x(i)), g(x(i)), q(x(i)), and g′(x(i)).The objective is given in equation (11a).Source: G Finite-dimensional Classical Stein… · p. 12
• Adjacent-point constraints bound derivatives by c2 and function increments by c2 times the corresponding point spacing.These restrictions appear in equation (11b) for every consecutive sorted pair.Source: G Finite-dimensional Classical Stein… · p. 12
• The notation includes positive-part and boundary-related quantities defined using c3.Equation notation specifies (r)+ ≜ max(r, 0) and defines Lb and Lu through c3.Source: G Finite-dimensional Classical Stein… · p. 12
• The optimization is finite-dimensional because vectors of function values and derivatives at sample or boundary points suffice.Introducing slack variables yields a convex quadratically constrained quadratic program with O(n) constraints, O(n) variables, and a linear objective.Source: G Finite-dimensional Classical Stein… · p. 12
• Extensions preserve the Stein objective because it depends on g and g′ only through their values at sample points.This observation underpins the equivalence between the finite program and the function-class formulation.Source: G Finite-dimensional Classical Stein… · p. 12

• A c2-bounded, c3-Lipschitz derivative must satisfy the program’s sampled derivative constraints.The proof considers such extensions and checks that equation (11b) follows.Source: §G.1 · p. 12
• The lower and upper derivative envelopes B and U bound every feasible derivative throughout each interval.Violating either envelope would contradict the c3-Lipschitz property relative to a sample point.Source: §G.1 · p. 12
• Integrating the derivative envelopes produces the lower and upper function-increment constraints in equations (11c) and (11d).The proof identifies the resulting integral bounds with the program’s constraints.Source: §G.1 · p. 12
• Boundary points satisfy the endpoint condition because the function is zero on the boundary.For interior points, the proof assumes a nonnegative derivative without loss of generality and integrates around the point.Source: §G.1 · pp. 12–13
• Integrating the Lipschitz derivative lower bound yields the quadratic constraint linking g(x(i)) and g′(x(i)) to c1.The argument obtains 1/(2c3) g′(x(i))^2 + g(x(i)) ≤ c1, with an analogous bound on the other side.Source: §G.1 · p. 13

G.2 Extending Feasible Solutions

• Every feasible program solution is extended by constructing lower and upper derivative functions independently on each adjacent interval.The construction requires c2-bounded, c3-Lipschitz functions matching g′ at both interval endpoints.Source: §G.2 · p. 13
• The extension derivative is a convex combination of the interval-wise bounds, with a coefficient chosen to match the required function increment.The construction defines the extension from a pair (mi, Mi) and a ζi ∈ [0, 1].Source: §G.2 · p. 13
• Convexity preserves c2-boundedness and c3-Lipschitzness, so only function boundedness at endpoints and critical points requires verification.The program’s constraint (11e) bounds endpoints, reducing the remaining check to points where the extension derivative vanishes.Source: §G.2 · p. 13
• The envelope functions B and U match endpoint derivatives and bound the derivative throughout each interval, reducing the construction to three sign cases.They serve as the building blocks for selecting mi and Mi.Source: §G.2 · p. 13
• When B and U retain one sign, choosing mi = B and Mi = U makes the extension monotonic and therefore c1-bounded.Monotonicity ensures boundedness for any convex-combination coefficient ζi.Source: §G.2 · p. 13
• When envelopes change signs, roots partition the interval and case-specific choices of mi and Mi control critical-point values within c1.The construction handles one-sign changes using φ and both-sign changes using ψ, with root-based bounds.Source: §G.2 · pp. 13–14
• The two-sign-change construction selects between B, W, and U according to the increment, preserving boundedness or well-definedness at the relevant root.For increments at most ψ, the construction uses B and W; otherwise it uses W and U.Source: §G.2 · p. 14
• The same construction extends the function beyond finite sample extremes when α or β is infinite, while retaining c2-bounded, c3-Lipschitz derivatives and |g| ≤ c1.The proof explicitly treats the left tail and states that the right tail is analogous.Source: §G.2 · p. 14

• The equivalence applies when X is a rectangular support and G_t is a t-spanner of supp(Q).Source: H Equivalence of Constrained Classical… · p. 14
• A t-spanner graph yields a graph Stein discrepancy bounded below by the classical discrepancy.Source: H Equivalence of Constrained Classical… · p. 14
• The graph Stein discrepancy is at most t^2κ_d times the classical Stein discrepancy.Here, κ_d depends only on dimension d and is independent of Q, P, G_t, and t.Source: H Equivalence of Constrained Classical… · p. 14

• The first inequality follows from containing the classical Stein set within the constrained graph Stein set.Source: Proof · pp. 14–15
• For the reverse direction, each constrained test function is extended to a classical Stein function while preserving its values and gradients on supp(Q).The construction enforces zero boundary values and zero off-coordinate derivatives on the coordinate boundary sets.Source: Proof · p. 15
• The Whitney-Glaeser extension theorem produces a classical Stein function with dimension-dependent factor t^2κ_d and matching Stein discrepancy objective value.Source: Proof · p. 15
• The extension defines boundary derivatives through a minimum involving the Lipschitz penalty ζ∥z−b∥_1.The resulting derivative remains ζ-Lipschitz on the support and boundary sets.Source: Proof · p. 15
• Triangle, Hölder, boundary-compatibility, and Lipschitz arguments establish the required regularity bounds before extension.Source: Proof · pp. 15–16