Stein Self-Repulsive Dynamics: Benefits from Past Samples
Authors: Mao Ye
Tongzheng Ren
Qiang Liu
Venue:  (2020)
Keywords: Approximate Inference
Markov Chain Monte Carlo
Stein Variational Gradient Descent

Abstract: We propose a new Stein self-repulsive dynamics for obtaining diversified samples from intractable un-normalized distributions. Our idea is to introduce Stein variational gradient as a repulsive force to push the samples of Langevin dynamics 
away from the past trajectories. This simple idea allows us to significantly decrease the auto-correlation in Langevin dynamics and hence increase the effective sample size. Importantly, as we establish in our theoretical analysis, the asymptotic stationary distribution remains correct even with the addition of the repulsive force, thanks to the special properties of the Stein variational gradient. We perform extensive empirical studies of our new algorithm, showing that our method yields much higher sample efficiency and better uncertainty estimation than vanilla Langevin dynamics.

Under review as a conference paper at ICLR 2020
STEIN SELF-REPULSIVE DYNAMICS: BENEFITS FROM PAST SAMPLES
Anonymous authors Paper under double-blind review
ABSTRACT
We propose new Stein self-repulsive dynamics for obtaining diversified sam- ples from intractable un-normalized distributions. Our idea to introduce Stein variational gradient ics away from the past trajectories. This simple idea allows us to significantly decrease the auto-correlation in Langevin dynamics and hence increase the ef- fective sample size. Importantly, as we establish in our theoretical analysis, the asymptotic stationary distribution remains correct even with the addition of the repulsive force, thanks to the special properties of the Stein variational gradient. We perform extensive empirical studies of our new algorithm, showing that our method yields much higher sample efficiency and better uncertainty estimation than vanilla Langevin dynamics.
repulsive force to push the samples of Langevin dynam-
INTRODUCTION
Drawing samples from complex un-normalized distributions one of the most basic problems in statistics and machine learning, with broad applications to enormous research fields that rely on probabilistic modeling. Over the past decades, large amounts of methods have been proposed for approximate sampling. including both Markov Chain Monte Carlo (MCMC) (e.g., Brookss et al.) 2011) and variational inference (e.g., Wainwright et al. 2008).
MCMC works by simulating Markov chains whose stationary distributions match the distributions of interest. Despite nice asymptotic theoretical properties, MCMC is widely eriticized for its slow convergence rate in practice. In difficult problems, the samples drawn from MCMC are often found to have high auto-correlation across time, meaning that the Markov chains explore very slowly in the configuration space. When this happens, the samples returned by MCMC only approximate small local region, and under-estimate the probability the regions un-explored by the chain.
Stein variational gradient descent (SVGD) Liu & Wangl 2016) is another type of approximation sampling methods designed to overcome the limitation of MCMC. Instcad of drawing random sam- ples sequentially, SVGD evolves pre-defined number of particles (or sample points) in parallel with a special interacting particle system to match the distribution of interest by minimizing the KI. divergence. In SVGD, the particles interact with each other to simultaneously move towards the high probability regions following the gradient direction, and also move away from each other due special repulsive force. As result, SVGD allows represent the variation of the distribution of interest. SVGD has been found promising tool for solving difficult sampling problems which diversity promotion critical (c.g../engg et al. 2017
obtain diversified samples that correctly
Haarnoja et al. 2017 Pu et al. 2017 Liu et al. 2017 Gong et al. 2019)) Various extensions have been developed (e.g.. Han & Liu 2018 Chen et al. 2018 Liu et 2019 Wang al. 2019a).
However, one problem of SVGD is that it theoretically requires to run infinite number of chains in parallel in order to approximate the target distribution asymptotically (Liu) 2017). With a finite number of particles, the fixed point of SVGD does still provide prioritized, partial approximation to the distribution in terms of the expectation Nevertheless. it is still desirable to develop variant of "single-chain SVGD". which only requires to run a single chain sequentially like MCMC to achieve the correct stationary distribution asymp- totically in time, with need to take the limit of infinite number of parallel particles.
a special case of functions (Liu Wang 2018).

---

Under review as a conference paper at ICLR 2020
In this work, we propose example of single-chain SVGD by integrating the special repulsive mechanism of SVGD with gradient- based MCMC such as Langevin dynamics. Our idea to use repulsive term of SVGD to enforce the samples in MCMC away from the past samples visited at previous iterations. Such a new self-repulsive dynamics allows us to decrease the auto-correlation in MCMC and hence increase the mixing rate, but still obtain the same stationary distribution thanks to the special property of the SVGD repulsive mechanism.
We provide throughout theoretically analysis of our new method, establishing it asymptotic conver- gence to the target distribution. As we show in the work, the analysis is highly non-trivial, because our new self-repulsive dynamics is a non- non-linear, high-order Markov process. Empirically, we ex- tensively evaluate our methods on array of challenging sampling tasks, showing that our method yields much better uncertainty estimation and larger effective sample size.
BACKGROUND: LANGEVIN DYNAMICS AND SVGD
In this section, we give brief introduction on Langevin dynamics (Rossky et al. 1978 and Stein Variational Gradient Descent (SVGD) (Liu & Wang 2016), which we integrate teogehter to develop our new self-repulsive dynamics for more efficient sampling.
Langevin dynamies Langevin dynamics is basic gradient based MCMC method. For some target distribution on Rd with density function (0) x exp(-V(0)), where V: Rd + R is the potential function, the (Fuler-discrerizee) Langevin dynamics simulates Markov chain with the following rule:
0x+1 = 0 (VV((e) + 2nex, CN,,)
where & denotes the number of iterations, }ex} are independent standard Gaussian noise, and rj is a step size parameter. It is well known that the limiting distribution of Og. when k oo approximates the target distribution when is sufficiently small.
Because the updates in Langevin dynamics local and incremental, points generated by the dynamics is highly correlated to the past sample. As result, we need to run Langevin dynamics sufficiently long order obtain diverse samples.
Stein Variatinal Gradient Descent (SVGD) Different from Langevin dynamics, SVGD itera- tively evolves a pre-defined number of particles parallel. Starting from an initial set of particles (0): i i= I, M) SVGD updates the M particles in parallel by
0 0 0 + ( 8M) Vi === M,
where the veloeity field, which denote by g( ). depends the empirical distribution of the current set of particles 8M == 314 Em oøg the following way.
g00: 8N = B
K(0, ((0.0))VV( + VoK O
Confining Term
Repulsive Term
Here dg is the Dirac measure centered at 0, and hence Eg [|] denotes averaging on the particles. The K(.,)) )) is positive definite kernel, such as the Gaussian RBF kernel, specified by users.
Note that g(0 ) consists of confining term and repulsive term: the confining particles to move towards high density region, and the repulsive term prevents the particles from colliding with each other. It is the balance of these two terms that allows us to asymptotically approximate the target distribution p²(0) x exp(-V(()) at the fixed point, when the number of particles goes to infinite. We refer the readers Liu & Wang (206)) Liu (011)); Liu Wang (2018) for throughout theoretical justifications of SVGD. But quick, informal way to justify the SVGD update is through the Stein identity, which shows that the velocity field g(;; p) equals zero when p equals the true distribution p", that is.
pushes
g(0'; :p") Ea-p' K((0 V'VV VoK(B,0)) =0, AO' Rd.
2

---

Under review as a conference paper at ICLR 2020
This shows that SVGD would converge if the particle distribution already forms closed approxi- mation to the target distribution p. meaning that the target distributions forms an (approximate) fixed point of the update.
STEIN SELF- REPULSIE DYNAMICS
In this work, we propose to integrate Langevin dynamics and SVGD to simultaneously decrease the auto-correlation of Langevin dynamics and climinate the need for running parallel chains in SVGD. The idea is to use Stein repulsive force between the the current sample and the past samples, hence forming a new self-avoiding dynamics with fast convergence speed.
Specifically. assume particle at the k-tth iteration. Denote by 8, the empirical distribution of M samples taken from the past iterations, i.e.,
single Markov chain like Langevin dynamics, where 0, denotes the
M Σ
60
Gn = =/11/
M
j=1
where Cy a thinning factor, which scales inversely with the step size 7) introduced to slim the sequence of past samples. Compared with the 8M in SVGD, which is averaged over M parallel particles, 8dd is averaged across time over M past samples. Given this, our Stein self-repulsive dynamics updates the sample via
0 + 0s + (-n)((0)) + 2nee.) + nag(Og; JM
Langevin
Stein Repulsive
in which the particle is updated with the typical Langevin gradient, plus Stein repulsive force against the samples from the previous iterations. a 2 0 is parameter that controls the magnitude of the Stein repulsive term. In this way, the particles pushed away from the past samples, and hence admits lower auto-correlation and faster convergence speed. Importantly, the addition of the repulsive force does not impact the asymptotic stationary distribution, thanks to Stein's 'ii's identity in (I). This is because when if the self-repulsive dynamics have converged to the target distribution p", such that 6₂ ~ p" for all k, the Stein self-repulsive term would equal to zero in expectation due to Stein's 's identity and hence does not introduce additional bias over Langevin dynamics. Rigorous theoretical analysis of this idea is developed in Section
Practical Algorithm Because B.M is averaged across the past samples, it is necessary to introduce burn-in in phase with the repulsive dynamics. Therefore, overall procedure works follows,
O₂ (Ox) + /2
k < Mc,
0 a
O2 + - -V (0s) + og(0,; )/0) + /2TEE
k 2 Mcg-
includes two phases. The first phase is the same as the Langevin dynamics which collects the initial M samples used in the second phase while serves warm start. The repulsive gradient update is introduced in the second phase to encourage the dynamics to visit the under-explored density region. We call this particular instance of our algorithm Self- Self-Repulsive Langevin dynamics (SRLD), self-repulsive variants of more general dynamics is discussed in Section 5
Remark Notice that, the first phase is introduced collect the initial M samples. However, it's not really necessary to generate the initial M samples with Langevin dynamics. We can simply use some other initialization distribution and get M initial samples from that distribution. In practice, we find using Langevin dynamics to collect the initial samples is natural and it can also be viewed as the burn-in phase before sampling, so we use all of the other experiments.
Remark The general idea of introducing self-repulsive terms inside MCMC other iterative algorithms is not new itself. For example, in molecular dynamics simulations, an algorithm called metadynamics (Laio & Parrinello, 2002) has been widely used. in which the particles are repelled away from the past samples in a way similar to our method, but with a typical repulsive function, such Σ; D) (Ox; Ok-je,,): where D(.,)) is any notation of dis-similarity. However, introducing
3

---

Under review as a conference paper at ICLR 2020
Hut
Lasgevis inttall aample
desalty
sumple
Combiaed Lomitty
Initial Sample
Initial Sample
Repulsive Sample
Langevin Sample
Figure I: An illustrative example showing the advantage of our Self-Repulsive Langevin dynamics.
With set of initial examples locating on the left part of the target distribution (show in yellow), Self-Repulsive Langevin dynamics is forced explore the right part more frequently, yielding an accurate approximation when combined with the initial samples. Langevin dynamics. however, does not take the past samples into account and yields poor overall approximation.
an arbitrary repulsive force would alter the stationary distribution of the dynamics, introducing a harmful bias into the algorithm. The key highlight of our approach, as reflected by theoretical results in Section [ is the unique property of the Stein repulsive term, that allows us to obtain the correct stationary distribution even with the addition of the repulsive force.
Remark Recently, (Zhang et al. 2018) proposed different combination of SVGD and Langevin dynamics, in which the Langevin force is directly added to set of particles that evolve in parallel with SVGD. Using our terminology, their system
0 0 (-))((O)) 22ee)) nog((); BM), ex IO
=i= I,..,
This is significantly different from our method on both motivation and practical algorithm. Their algorithm still requires to simulate M chains of particles in parallel like SVGD, and was proposed to obtain easier theoretical analysis than the deterministic dynamics of SVGD, Our work instead motivated by the practical need of decreasing the auto-correlation in Langevin dynamics, and avoid- ing the need of running multiple chains SVGD, and hence must be based on self-repulsion against past samples along single chain.
In another recent work, (Chen et al. 2018 proposed a ग T-SLD method, which simulates the lin- combination of the evalutionary partial differential equations of SVGD and Langevin dynamics using discrete gradient flow with blob-hased method. Their method again motivated by the theo- retical interest of discovering categories of algorithms, and does not involve self-repulsive single chain like our method.
An Illustrative Example Here we give an illustrative example to show the key advantage of our self-repulsive dynamics. Assume that we want to sample from a bi-variate Normal distribution shown in Figure [] Unlike standard settings, we assume that we have already obtained some initial samples (yellow dots in Figure T) before running the dynamics. The initial samples are assumed to concentrate on the left part of the target distribution shown in Figure [ In this extreme case, since the left part the distribution is over-explored by the initial samples, it is desirable to have the subsequent new samples to concentrate more the un-explored part of the distribution. How- ever, standard Langevin dynamics does not take this into account, and hence yielding a bias overall representation of the true distribution (see the left panel). With self self-repulsive dynamics, the new
4

---

Under review as a conference paper at ICLR 2020
Practical Dynamics
Theorem 4.3
Discrete- Time
Theorem 4.2
Continuous- Time Mean-Field Limit
Theorem 4.1
Stationary
Mean-Field Limit
Distribution
ttPartclles
Step size 11 -
p
M 00
Figure 2: Roadmap of the thcoretical analysis. Theorem 4.3 shows the mean-field limit when M + oo.
Theorem 4.2 bounds the time discretization error. And Theorem [.T shows that the limiting distribution of the continuous-ime mean field dynamics the target distribution want.
samples are forced to explore the un-explored region more frequently, allowing us to obtain much more accurate approximation when combining the new and initial samples.
THEORETICAL ANALYSIS OF STEIN SELF-REPULSIVE DYNAMICS
In this section, we provide theoretical analysis of the self-repulsive dynamics. We establish that our self-repulsive dynamics converges to the correct target distribution asymptotically, in the limit when M approaches to infinite and the step size 7 approaches to O. This is highly non-trivial task, because the self-repulsive dynamics is a highly complex, non-linear and high order Markov stochastic process. We attack this problem by breaking the proof into the following three steps
illustrated in Figure 2
1) At the limit of large particle sizes M + (called the mean field limit), show that practical dynamics in (3) is closely approximated by a discrete-time mean-field dynamics characterized by below.
2) By further taking the limit of small step size + 0+ (called the continuous-time limit), the dynamics in converges a continuous-time mean-field dynamics characterized by below.
3) We show that the dynamics in S converges to the target distribution.
Remark As mentioned in Section 3 introduce the first phase collect the initial M sam- ples for the second phase, and our theoretical analysis follows this setting to make our theory as close the practice. However, the theoretical analysis can be generalized to the setting of drawing M initial samples from some initialization distribution with almost identical argument.
Notations We use |-|l and (:-)) represent the & vector norm and inner product respectively. The Lipschitz norm and bounded Lipschitz norm of a function f are defined by Il and [[] The KL divergence, Wasserstein-2 distance and Bounded Lipschitz distance between distribution Pl, P2 are denoted DKL/P1 P2], W2 P2] and [BL[PI Pa] respectively.
4.1 MEAN-FIELD AND CONTINUOUS TIME LIMITS
To fix the notation, we denote by == Law(Or 0g) the distribution of Og at time i of the practical self-repulsive dynamics (3), which we
we refer as the practical dynamics in the sequel, when the initial particle 0, is drawn from an initial continuous distribution Po supported on Rd. Note that given Po, the subsequent Ak can be recursively defined through dynamics (3). Due the diffusion noise in
Langevin dynamics, all P are continuous distributions supported on I We now introduce the limit dynamics when we take the mean-field limit (M + +o0) and then the continuous time limit (n + 0+).
Discrete- Rismrete-Time Mean-Field Dynamies (M + +0o) In limit M - oo, we show that our practical dynamics is closely approximated by the following dynamics, in which the empirical measures are replaced by the continuous distributions the particles themselves.
0 vV.) + /2TeR.
t < Mc
0.+1 =
(4)
0 + (V)))) + (Ak. mH1) + /2nex,
1 2 Mc.
where AM = M Zj=1 Pk-jan by A = Law Law(Oz) the (smooth) distribution of 0 at time-step k when the dynamics is initialized with O go Pa- Compared with the practical dynamics
M
5

---

Under review as a conference paper at ICLR 2020
in (3). the difference is that the empirical distribution 8 is replaced by the smooth distribution 7 Similar to the recursive definition of Pk. following dynamics (3) Pb. is also recursively defined through dynamics ( starting from Lo = Po- As we show Theorem 4.3) if the auto-correlation of O2 decays fast enough and M is sufficiently large, p is well approximated by the empirical
distribution MM in (3). and further the two dynamics ((3) and (4) converges to each other in the sense that W2lPk Pk] + 0 as M - oo for any k. Note that in taking the limit of M, we need to ensure that we run the dynamics for more than Mc, steps. Otherwise, SRLD degenerates to Langevin dynamics as we stop the chain before we collect M samples.
Continuous- Time Mean-Field Dynamics (1) - 0+) In the limit of zero step size (1) + 0+), the discrete-time mean field dynamies in (4) can be shown to converge to the following continuous-time mean-field dynamics:
d0, VV(0,)dt + dB,, VV(O.) ag(Ok, ) dt dB,, 2 Mc,
t E [0,
(5)
where P. M ΣØ Pt=10t-j(()) and Pt = Law (0.) is the distribution of 0, at a continuous time point with G, initialized by 0g गU Po. We prove that (5 is closely approximated by (4) with small step size in the sense that [Px. || pknpl + 0 as + O Theorem 4.2 and importantly, the stationary distribution of (5) equals to the target distribution p²()) x exp(
(0))
4.2 ASSUMPTIONS
We first introduce the techinical assumptions used in theoretical analysis. Assumption 4.1. (RBF Kernel)
use RBF kernel ie.. K(O,, 02) = exp(- |00 - 02112 /o) for some fixed a < oo.
We only assume the RBF kernel for the simplicity of our analysis. However, is straightforward to generalize theory to other positive definite kernels.
Assumption 4.2. (V dissipative and smooth)
Assume that (0. (V) S a1 |101 and - s b1 110 02ll. We also assume that VV]] S by. Here a1 and b are some finte positive constant.
Assumption 4.3. (Regularity Condition)
Assume [1 < O. Define pM ΣØ= Pk- Neejje./M M, assume there exists B < such that
9((0 85/ g(0; pb 1²)
inf 12 M2
sup E
|
Assumption 4.4. (Strong-convexitty
Suppose that -VV(O))-V(00)).)) 02) > L L| - -211 for positive constant L.
Remark Assumption is standard in the existing Langevin dynamics analysis (see Dalalyan Raginsky et al. (211)) for example). Assumption 4.3 is weak condition as it assumes that the dynamics can not degenerate into one local mode and stop moving anymore. This assumption expected to be true when we have diffusion terms like the Gaussian noises in our self-repulsive dynamics. Assumption 4.4 is classical assumption on the existing Langevin dynamics analysis with convex potential [alalyan] (0017); Durmus et al. (2019). Although being a bit strong, this assumption broadly applies to posterior inference problem in the limit of big data, as the posterior distribution converges to Gaussian distributions for large training set as shown by Bernstein -von von Mises theorem. It is technically possible to further generalize results to the non-convex settings with a refined analysis, which we leave as future work. This work focuses on the classic convex setting for simplicity.

---

Under review as a conference paper at ICLR 2020
4.3 MAIN THEOREMS
All of the proof in this section can be found in Appendix B.5
We first prove that the limiting distribution of the continuous-time mean field dynamies 5 target distribution. This is achieved by writing dynamics into the following (non-linear) partial differential equation:
is the
V (-VVp)) A.p. depr V [(-VV + (. pc] Spc. 2 Me.
t [0, Mc)
Theorem 4.1. (Stationary Distribution)
Given some finite M, e and o, and suppose that the limiting distribution of dynamics 5
exists. Then
the limit distribution is unique, and equals to (0) x exp(-V(0)).
We then give the upper bound on the discretization error, which can be characterized by analyzing the KI. divergence between PE and Pken
Theorem 4.2. (Time Discretization Error)
Given some sufficiently small step size n and choose a < 02/( + 4/0). Under assumption 4.3 and Cy c/n/ we have for some constant C,
10 (n (0(++hi)) knf² {0
k < Mcn -1
DKL [pin P1] 1E(O....k)
n Meq + a² a Oa
k Mep-
With this theorem, we can know that if is small enough, then the discretization error is small and I approximates p closely. Next we give result on the mean field limit of M 00.
Theorem 4.3. (Mean-Fiee) Field Limit)
4.2 and suppose that we choose a and n such that (a1 2ody /o) mb; 20 (b, + +1
Under assumption 4.1
< 1; 02 a (2b) =) > 0: Then there exists a constant C2. such that when
L/a 2 02 and we have
Welp Pk = o 0 (a²/M+ for k Mcn.
and W ( 0ff < Mcy
Remark Combine all the result, we have lim Dar [e,, p-] = O, where the joint limit requires kn exp(Caa y S 1, the dynamics is degenerated to Langevin dynamics.
o(1) and Me
kn
(1 o(1)), with y > 1. Notice that if
EXTENSION TO GENERAL DYNAMICS
Although we have focused on self-repulsive Langevin dynamics, our Stein self-repulsive idea can be broadly combined with general gradient-based MCMC. Following Ma et al. (2015), we consider the following general class of sampling dynamics for drawing samples from p(0) x exp (exp(-V/0 (())):
d0, f(0)dt + (2D(O),B,,
0 (Di, (0) (0). Σ j=1 00
with =(0)= (D(0) Q())VV r(0),
r.(0)
where D is positive semi-definite diffusion matrix that determines the strength of the Brownian motion and Q is a skew-symmetric curl matrix that can represent the traverse effect (c.g. in NNal al. 2011 Ding et al. 2014). By adding the Stein repulsive force, obtain the following general self-repulsive dynamics
(- f(0)dt dB, (f(0) + (g; PM dt + dB, t Me
E [0, Mc)
d0, {
(6)
7

---

Under review as a conference paper at ICLR 2020
MMD
Wasserstein Distance
Effective Sample Size
Autocorrelation
0.12
Repuisive Dynamics
25
1.0
Langavin Dynamies
0.08
0.8
20
15
0.4
0.04
10
0.1
0.00
0.2
100
780
1300
1900
100
700
1300
1900
100
700
1300
1900
200
400
600
505
Num of samples
Num of samples
Num samples
Num of lags
Figure 3: Sample quality of SRLD and Langevin dynamics for sampling the correlated 2D distribu- tion. The auto-correlation the averaged auto-correlation of the two dimensions.
where P: = Law (e:) is again the distribution 0, following O when initalized at 0g ~ Pg. Similar to the case of Langevin dynamics, this process also converges to the correct target distribution, and can be simulated by practical dynamics similar to Theorem 5.1. (Stationary Distribution)
Given some finite M, e and o,, and suppose that the limiting distribution of dynamics (6) exists. Then the limiting distribution is unique and equals the target distribution p" (0) exp( expp-V((0) (())..
6 EXPERIMENTS
In this section, we evaluate the proposed method various challenging problems, including sam- pling the posteriors of Bayesian Neural Networks, and uncertainty estimation in Reinforcement Learning. Our results show that our Self-Repulsive Langevin dynamics (SRLD) yields much higher sample efficiency than vanilla Langevin dynamics. Our code is available along the submission.
6.1 SYNTHETIC EXPERIMENT
This experiment aims to show how the repulsive gradient helps explore the whole distribution. Following Ma al. (2015), we compare the sampling efficiency on the following correlated 2D distribution with density
A Correlated 2D Distribution
([01, 02)) / (02 + - 03)
We compare the SRLD with vanilla Langevin dynamics, and evaluate the sample quality by Maxi- mum Mean Discrepancy (MMD) (Gretton et al. 2012), Wasserstein- Distance and effective sample size (ESS). Notice that the finite sample quality of gradient based MCMC method is highly related to the step size. Compared with Langevin dynamics, we have an extra repulsive gradient and thus implicitly have larger step size. To rule this effect, we different step sizes of the two dy- namics so that the gradient of the two dynamics has the same magnitude. In addition, to decrease the influence of random noise, the two dynamics are set to have the same initialization and use the same sequence of Gaussian noise. We collect the sample of every iteration. We repeat the experiment 20 times with different initialization and sequence of Gaussian noise.
Figure 3 summarizes the result with different metrics. We can see that SRLD have signifi- cantly smaller MMD and Wasserstein-1 Distance as well as larger ESS compared with the vanilla Langevin dynamics. Moreover, the introduced repulsive gradient creates negative auto-correlation between samples. sigure4ssows a typical trajectory of the two sampling dynamics. We can see that SRLD have faster mixing rate than vanilla Langevin dynamics. Note that since sequence Gaussian noise for both algorithms, the difference is mainly due to the use of repulsive gradient rather than the randomness.
the same
Mixture of Gaussian Distribution We aim to show how the repulsive gradient helps the particle escape from the local high density region by sampling the 2D mixture of Gaussian distribution using SRLD and Langevin dynamics. The target density is set to be
1 p" p()))= x 2 exp (-11- - 112 +
0 + 11/ 12), p= =10,, 101,0z] 0z]T.
exp
8

---

Under review as a conference paper at ICLR 2020
Repulsive Dynamics
Langevin Dynamics
Num of samples 50
Nom of samples 100
Num of samples = 180
Num samples 200
Num of sumples 250
Figure 4: Sampling trajectory of the correlated 2D distribution.
MMD
Wasserstein Distance
Effective Sample Size
A Autocorrelation
28
0.16
Repulsive Dynamics
Langovin Dynamies
1.93
24
0.12
20
16
0.08
1.09
0.04
0.67
0.00
100
300
500
700
900
100
300
900
100
300
500
700 Num of lags
000
800
Num of samples
Num of samples
Num of samples
Figure 5: Sample quality and autocorrelation of the mixture distribution. The auto-correlation is the averaged auto-correlation of the two dimensions.
where [1, This target distribution have mode at -1 and 1, and vanilla Langevin dynamics can stuck in one mode while keeps the another mode under-explored (as the gradient of energy function can dominate the update of samples). We use the same evaluation method, step sizes, initialization and Gaussian noise as the previous experiment. We collect one sample every 100 iterations and the experiment repeated for 20 times. HigureB]show that SRLD again consistently outperforms the Langevin dynamics on all of the evaluation metrics.
6.2 BAYESIAN NEURAL NETWORK
Bayesian Neural Network is one of the most important methods in Bayesian Deep Learning with wide application in practice. Here we test the performance of SRLD on sampling the posterior of Bayesian Neural Network on the UCI datasets (Dua & Graff 2017). We assume the output normal distributed, with a two-layer neural network with 50 hidden units and tanh activation to predict the mean of outputs. We set [[1, 0.1) prior for the inverse output variance. All of the datasets randomly partitioned into 90% for training and 10% for testing. The results are averaged over 20 random trials. We the mini-batch size be 100 and the number of past samples M to be 10. In all experiments, we use RBF kernel with bandwidth set by the median trick as suggested in Liu Wang (016)) For SVGD, (2016). We run 50000 iterations for each methods, and for LD and SRLD, the first 40000 iteration discarded as burn-in. We use a thinning factor of C, c/n 100 and in total we collect 100 samples from the posterior distribution. For each dataset, we generate 3 extra data splits for tuning the step size for each method.
the original implementation with 20 particles by Liu & Wang
Table [ shows the average test RMSE and test log-likelihood and their standard deviation. The method that has the best average performance marked as boldface. We observe that large portion of the variance is due to the random partition of the dataset. Therefore, show the statistical significance, we use the matched pair t-test to test the statistical significance, mark the methods that perform the best with significance level of 0.05 with underlines. Note that the results of SRLD/LD and SVGD is not very comparable, because SRLD/LD are single chain methods which averages across time, and SVGD multi-chain method that only the results of the last iteration. We

---

Under review as a conference paper at ICLR 2020
Ave Test RMSE
Ave Test L.
Dataset
SVGD
LD
SRLD
SVGD
LD
SRLD
Boston
3.300 0.111
3.312
0.187
0.188
-1.276
-2,778
-2,500 0.054
Concrete
4.994 1: 0.171
4.908
0.113
4.886 1 0.108
-5.500 0.399
-3,055 4 0.035
-3.034
0.031
Energy
0.428 = 0.016
0,412
= 0.016
0.305
-0.781
-0.543 0.011
-0.476+0.030
Naval
0,006 +: 0.000
0.006
0.002
1.00000.000
3.056
0.034
= 0,030
4.186
0.015
WineRed
0.655 + 0.008
0.649 = 0.009
0.639
+ 0.0
-1001
-1.001
0.019
-0.970
0.066
WineWhite
0.655
0.692
-1.040
-1.033
0.000
Yacht
0.593 + 0.071
0.597 0.051 0.578 + 0.054
-1.281 0.279
-1.187 # 0.307
0.458 0.036
Table Averaged test RMSE and test log- likelihoodd on UCI datasets. Results are averaged over 20 trails. The boldface indicates the method has the best average performance and the underline marks the methods that perform the best with a significance level of 0.05.
provide additional results in appendix that SRLD averaged on 20 particles (across time) can also achieve similar or better results as SVGD with 20 (parallel) particles.
6.3 CONTEXTUAL BANDITS
We evaluate the quality of uncertainty estimation provided by our methods on several contextual bandits problems. Uncertainty estimation is a key component of contextual bandits. If the agent makes decisions with a poorly estimated uncertainty, the decisions will finally lead catastrophic failure through the feedback loops (Riquelme et al. 2018).
Though in principle all of the MCMC methods return the samples follow the true posterior if we can run infinite MCMC steps, in practice we can only obtain finite samples as we only have finite time run the MCMC sampler. In this case, the auto-correlation issue can lead to the ender-estimate the uncertainty, which will cause the failure on all of the reinforcement learning problems that need exploration.
We consider the posterior sampling (a.k.a Thompson sampling) algorithm with Bayesian neural network as the function approximator. We follows the experimental setting from Riiqelmme et al. (018). The only difference is that change the optimization of the objective (e.g. evidence lower bound (ELBO) in variational inference methods) into running MCMC samplers. We set the step of samplers equal to the optimization step, and use thinning factor of 100. We compare the SRLD with the Langevin dynamics on the Mushroom and Wheel bandits from (Riquelme et al. 2018). and include SVGD baseline. For more detailed introduction and setup, see Appendix AAS A.5
The cumulative regret is shown in Table 2 SVGD is known to have the under-estimated uncertainty for Bayesian neural network if particle number is limited (Wang et al. 2019b), and result, has the worst performance among the three methods. SRLD is slightly better than vanilla Langevin dynamics on the simple Mushroom bandits. On the much more harder Wheel bandits, SRLD is significantly better than the vanilla Langevin dynamics, which shows the improving uncertainty estimation of our methods within finite number of samples.
Dataset
SVGD
LD 4.28
SRLD 3.80
Mushroom
20.7 2.0
0.09
0.16
Wheel
91.32 0.17
38.07 + 1.11
32.08 + 0.75
Table 2: Cumulative Regrets on two bandits problem. Results are averaged over 10 trails. Boldface indicates the methods with best performance and underline marks the best significant methods with significant level 0.05.
CONCLUSION
We propose a Stein self-repulsive dynamics which applies Stein variational gradient to push samples from MCMC dynamics away from its past trajectories. This allows us significantly decrease the cottocoorelation of MCMC, increasing the sample efficiency for better estimation. The advantages of our method are extensive studied both theoretical and empirical analysis in our work. In future work, we plan to investigate the combination of our Stein self-repulsive idea with more general MCMC procedures, and explore broader applications.
10

---

Under review as a conference paper at ICLR 2020
REFERENCES
HZ An and FC Huang. The geometrical ergodicity of nonlinear autoregressive models. Statistica Sinica, pp. 943 936,956, 1996.
Charles Blundell. Julien Cornebise, Koray Kavukcuoglu, and Daan Wierstra. Weight uncertainty in neural network. In Francis Bach and David Blei (eds.), Proceedings of the 32nd International Conference Machine Learning, volume 37 of Proceedings of Machine Learning Research, pp. 1613-1622, Lille, France, 07-09 Jul 2015. PMLR. URL htttt:::/pprroceedngg.. .ml.....s.........mmmmm........m
Richard Bradley. Basic properties of strong mixing conditions, a survey and some open questions. Probability Surveys, 2:107-144, 2005.
Richard C. Bradley. Introduction to strong mixing conditions. 2007.
Steve Brooks, Andrew Gelman, Galin Jones, and Xiao-Li Li Meng. Handbook of Markov chain Monte Carlo. CRC press, 2011.
Sébastien Bubeck, Nicolo Cesa- Besaa-Bincch,, et al. Regret analysis of stochastic and nonstochastic multi- armed bandit problems. Foundations and Trends(R) in Machine Learning. 5(1):1-122, 2012.
Olivier Chapelle and Lihong Li. An empirical evaluation of Thompson sampling. In Advances in neural information processing systems, pp. 2249-2257, 2011.
Changyou Chen, Ruiyi Zhang, Wenlin Wang, Bai Li, and Liqun Chen. A unified particle- optimization framework for scalable Bayesian sampling. In Proceedings of the Shirty-Fourdh Conference on Uncertainty in Artificial Intelligence, UAI 2018, Monterey, California, USA, Au- gust 6- 6-10, 2018, pp. 746-755, 2018.
Arnak S Dalalyan. Theoretical guarantees for approximate sampling from smooth and log-concave densities. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 79(3): 651-676, 2017.
Nan Ding, Youhan Fang, Ryan Babbush, Changyou Chen, Robert D. Skeel, and Hartmut Neven. Bayesian sampling using stochastic gradient thermostats. In Advances in Neural Information Processing Systems 27: Annual Conference on Neural Information Processing Systems 2014, December 8-13 2014, Montreal, Quebec, Canada, pp. 3203-3211, 2014.
Dheeru Dua and Casey Graff. UCI machine learning repository, 2017. Utt::tp//arcc ics.uci. uui.ucui.dum..
Alain Durmus, Szymon Majewski, and Blazej Miasojedow. Analysis of langevin monte carlo via convex optimization. Journal of Machine Leurning Research, 20(73):1-46, 2019.
Yihao Feng, Dilin Wang, and Qiang Liu. Learning draw samples with amortized stein variational gradient descent. uncertainty in artificial intelligence, 2017.
Chengyue Gong, Jian Peng, and Qiang Liu. Quantile Stein variational gradient descent for batch Bayesian optimization. In International Conference on Machine Learning, pp. 2347-2356, 2019.
Arthur Gretton, Karsten M Borgwardt, Malte Rasch, Bernhard Schölkopf, and Alexander Smola. A kernel two-sample test. Journal of Machine Learning Research, 13(Mar):723-773 2012.
Tuomas Haarnoja, Haoran Tang, Pieter Abbeel, and Sergey Levine. Reinforcement learning with deep ensrgy-based policies. international conference on machine learning, pp. 1352-1361, 2017.
Jun Han and Qiang Liu. Stein variational gradient descent without gradient. In Uncertainty in Artificial Intelligence, 2018.
Alessandro Laio and Michele Parrinello. Escaping free-energy minima. Proceedings of the National Academy of Sciences, 99(20): 12222)225622111 12562-12566666 2002.

---

Under review as a conference paper at ICLR 2020
Chang Liu, Jingwei Zhuo, Pengyu Cheng, Ruiyi Zhang, and Jun Zhu. Understanding and acceler- ating particle-base variational inference. International Conference on Machine Learning, pp- 4082-4092, 2019.
Qiang Liu, Stein variational gradient descent as gradient flow. In Advances in neural information processing systems, pp. 31 11-31226, 2017.
Qiang Liu and Dilin Wang. Stein variational gradient descent: A general purpose bayesian inference algorithm. Advances In Neural Information Processing Systems, pp. 2378-22-- 2386, 2016.
Qiang Liu and Dilin Wang. Stein variational gradient descent as moment matching, Advances in
Neural Information Processing Systems, pp. 8854 8863, 2018.
Yang Liu, Qiang Liu, and Jian Peng. Stein variational policy gradient. uncertainty in artificial intelligence, 2017.
Yi-An Ma, Tianqi Chen, and Emily Fox. A complete recipe for stochastic gradient memc. In Advances in Neural Information Processing Systems, pp. 2917-2925, 2015.
Radford M Neal et al. Meme using hamiltonian dynamics. Handbook of markov chain monte carlo, 2(11):2, 2011.
Yunchen Pu, Zhe Gan, Ricardo Henao, Chunyuan Li, Shaobo Han, and Lawrence Carin. Vae learn- ing via stein variational gradient descent. neural information processing systems, pp. 4236-4245, 2017.
Maxim Raginsky, Alexander Rakhlin, and Matus Telgarsky. Non-convex learning via stochastic gradient langevin dynamics: a nonasymptotic analysis. In Proceedings of the 30th Conference on Learning Theory, COLT 2017, Amsterdam, The Netherlands, 7-10 July 2017, pp. 1674-1703, html
2017. http::ttt////
mwr. wtt....p........//i.i.iaa....... hwtps//rgiiinkyy.y.
Carlos Riquelme, George Tucker, and Jasper Snoek. Deep bayesian bandits showdown: An empir- ical comparison of bayesian deep networks for thompson sampling. In 6th International Con- ference Learning Representations, ICLR 2018, Vancouver, BC, Canada, April 30 May 3, 2018, Conference Track Proceedings, 2018. URL https://openzevieeee net /for/ffrum id SyYe6k-KKS CW
Peter Rossky, JD Doll, and HI Friedman. Brownian dynamics as smart monte carlo simulation. The Journal of Chemical Physics, 69(10)46628 4633, 1978.
Jeff Schlimmer. Mushroom records drawn from the audubon society field guide to north american mushrooms. GH Lincoff (Pres), New York, 1981.
William R. Thompson. On the likelihood that one unknown probability exceeds another in view of the evidence of two samples. Biometrika, 25(3/4):255-944 1933.
Martin Wainwright, Michael I Jordan, et Graphical models, exponential families, and variational inference. Foundations and Trends in Machine Learning. 1(1-2):1-305, 2008.
Dilin Wang, Ziyang Tang, Chandrajit Bajaj, and Qiang Liu. Stein variational gradient descent with matrix-valued kernels. In Conference on Neural Information Processing Systems, 2019a.
Ziyu Wang, Tongzheng Ren, Jun Zhu, and Bo Zhang. Function space particle optimization for bayesian neural networks. arXiv preprint arXiv: 2019b.. 2019b.
Jianyi Zhang, Ruiyi Zhang, and Changyou Chen. Stochastic particle-optimization sampling and the non-asymptotic convergence theory. arXiv preprint wrXi::1809.01293,
12

---

E 0.5I. We run SRLD with 100, 50, 20, 10, 0 and the case collect sample every 10 iterations. The other experiment setting is the same as the toy examples set SRLD Under review as a conference paper at ICLR 2020
A EXPERIMENT DETAILS AND ADDITIONAL EXPERIMENT RESULT
A.I SYNTHETIC 2D MIXTURE OF GAUSSIAN EXPERIMENT
Republive Dynamics
Langevih Dynamics
Num of samples - 100
Num of samples 150
Num of samplas - 200
Num of samples 250
Num of samples - 300
Figure 6: Sampling trajectory of the mixture of Gaussian.
To provide more evidence on the effectiveness of SRLD escaping from local high density region, we plot the sampling trajectory of SRLD and vanilla Langevin dynamics on the mixture of Gaussian mentioned in Section 6.1 We can find that, when both of the methods obtain 200 samples, SRLD have started to explore the second mode. while vanilla Langevin dynamics still stuck in the original mode. When both of the methods have 250 examples, the vanilla Langevin dynamics just start to explore the second mode, while our SRLD have already obtained several samples from the second mode, which shows our methods effectiveness on escaping the local mode.
A.2 SYNTHETIC HIGIEE DIMENSIONAL GAUSSIAN EXPERIMENT
To show the performance of SRLD in higher dimensional case with different value of we addition- ally considering the problem on sampling from Gaussian distribution with d = 100 and covariance 0 reduces to Langevin. We
in the main text. The results summarized at Figure (7). In this experiment, with an inappropriate = 100. For this chain, the repulsive gradient gives strong repulsive force and thus has the largest ESS and the fastest decay of autocorrelation. While the inappropriate value induces too much extra approximation error and thus its performance is not as good as these with smaller a (see MMD and Wasserstein distance). This phenomenon matches our theoretical finding.
MMD
Wasserstein Distance
ESS
Autocorrelation
0.035
1.0
Alpba
20.0
120
Alpba- 20
Aigba 50
17.5
100
0.8
0.020
Alpha
80
0.015
A
12.5
60
0.010
10.0
7.5
0.0
100
300
500
700
900
100
300
500
700
900
100
300
500
700
900
12
15
Num of samples
Num of samples
Num of samples
Num of lags
Figure 7: Sample quality and autocorrelation of the higher dimensional Gaussian distribution. The auto-correlation is the averaged auto-correlation of all dimensions.
13

---

Under review as a conference paper at ICLR 2020
MMD
Wasserstein Distance
ESS
Autocorrelation
Alpla 20
12
0.015
Alpha
Appa-S
10
0.4
Alpba 0
0.3
15
0.2
0.005
10
0.1
0.000
500
700
100
300
100
300
900
20
Num of samples
Num of samples
Num of samples
Num of lags
Figure 8: Sample quality and autocorrelation of the higher dimensional mixture distribution. The auto-correlation is the averaged auto-correlation of all dimensions.
A.3 SYNTHETIC HIGHER DIMENSIONAL MIXTURE OF GAUSSIAN EXPERIMENT
We also consider sampling from the mixture of Gaussian with d 20. The target density is set to be
p²( x exp -0.5
2/d1 1) exp (-0.5 o v 2/d1
2
where 0 = (0.....000) T and 1 [I....,1] 1J And thus the mean of the two mixture component 20, 10, 5, (and when 0, it reduces to LD).
with distance 2/2. We run SRLD with The other experiment setting is the same as the low dimensional mixture Gaussian case. Figure (8] summarizes the result. As shown in the figure, when a becomes larger, the repulsive forces helps the sampler better explore the density region.
A.4 UCI DATASETS
We show some additional experiment result on posterior inference on UCI datasets. As mentioned in Section 6.2 the comparison between SVGD and SRLD is not direct as SVGD is a multiple-chain
method with fewer particles and SRLD is single chain method with more samples. To show more detailed comparison, we compare the SVGD with SRLD using the first 20. 40, 60, 80 and 100 samples, denoted as SRLD-n where r is the number samples used. Table 3 shows the result of averaged test RMSE and table 4 shows the result of averaged test loglikelihood. For SRLD with different number of samples, the value is set to be boldface it has better average performance than SVGD. If is statistical significant with significant level 0.05 using matched pair I-test, add an underline on it.
Figure 9 and 10 give some visualized result the comparison with Langevin dynamics and SRLD. To rule out the variance of different splitting on the dataset, the errorbar is calculated based on the difference between RMSE SRLD and RMSE of Langevin dynamcis in 20 repeats (And similarily for test liceliiod,.odd) And we only applied the error bar on Langevin dynamics.
Ave Test RMSE
Dataset
SRLD-20
SRLD-40
SRLD-60
SRLD-80
SRLD-100
SVGD
Beston
3.236 1 0.174
3.173 1 0.176
3.130 1 0.173
3.101
1. 0.179
3.086
0.181
3.300 1 0.142
Concrete
4.959
4.921
+ 0.111
4.906
0.109
4.891 0.108
4.886
00108
4.904
Energy
0.422 = 0.016
0.409
: 0.016
0.405
0.399 0.016
0.016
0.395
0.016
0.128
Naval
0.005 : 0.001
0.004 + 0.000
0.000.00.00
0.003
1: 0.000
0.003
0.000
0.006 - 0.000
WineRed
0.654 0.009
0.647
0.009
0.644
0.641
0.009
0.639
0.009
0.655 0.008
WineWhite
0.605 +: 0.003
0.692
1 0.003
0.690
1 0.003
0.689
1: 0.002
0.688
0.003
0.655 1 0.008
Yacht
= 0.055
0.008 + 0.052
0.597 0.051
0.587 0.054
0.578
0.054
0.593 0.071
Table 3: Comparing SRLD with different number of samples with SVGD on test RMSE. The results are computed over 20 trials. For SRLD, the value is set to be boldface if it has better average performance than SVGD. The value if with underline if it is significantly better than SVGD with significant level 0.05 using matched pair t-test.
14

---

Under review as a conference paper at ICLR 2020
Ave Test Test
SRLD-20
SRLD-40
SRLD-60
SRLD-80
SRLD-100
SVGD
-2.642 .088
-2.527
2.516100:88
-2.500
0.0051
-1.766 1.276
Concrete
-3.084 1: 0.038
-3.061 1. 0.034
-3,050
0.033
-3.040 1 0.031
-3.034
0.031
-5.500 1: 0.398
Energy
-0.580 + 0.033
-0.536 0.048
-0.522
0.016
-0,501+00011
-0.476
-0.781
0.094
Naval
4.033 1 0.230
4.100 +: 0.171
4.140
1.0011
4.167:0.011
3,056
WineRed
-1.008 + 0.019
-0.990 0.017
-0.982 = 0.016
-0.974
0.016
-0.970
0.016
-1,010 0.018
WineWhite
0.004
-1.043
-1.040 0.019
Yachs
1.160 = 0.2 -0,650 + 0.173
-0.556 0.096
-0.465 + 0.037
-0.458 + 0.036
-1.281 0.279
Table 4: Comparing SRLD with different number of samples with SVGD on test log-likelihood. The results are computed over 20 trials. For SRLD, the value is set to be boldface if it has better average performance than SVGD. The value if with underline it is significantly better than SVGD with significant level 0.05 using matched pair t-test.
Boston
Concrete
Energy
Naval
3.8
Repulsive Dynamics
5.15
0.70
0.035
Langevin Dynamics
0.030
3.6
0.65
0.025
5.05
0.60
3.4
5.00
0.55
0.015
3.3
0.50
4.95
0.010
3.1
4.90
0.40
0.005
20
40
60
60
100
0
20
40
60
80
100
20
40
60
80
100
0
20
40
60
80
100
Num of samples
Num samples
Num of samples
Num samples
WineRed
White
Yacht
0.75
0.69
0.74
0.725
0.700
0.67
0.650
0.66
0.600
0.64
0.575
0
20
40
60
60
100
20
40
60
60
100
0
20
40
60
80
100
Num samples
Num of samples
Num samples
Figure 9: Comparison between SRLD and Langevin dynamics on test RMSE, The results are com-
puted based on 20 repeats. The error bar is calculated based on RMSE of SRLD WALD-RMSE of Langevin dynamics in 20 repeats to rule out the variance of different data splitting
A.5 CONTEXTUAL BANDIT
Contextual bandit is class of online learning problems that can be viewed as a simple reinforcement learning problem without transition. For completely understanding of contextual bandit problems, we refer the readers to the Chapter 4 of (Bubeck al., (0122.. Here we include the main idea for completeness. In contextual bandit problems, the agent needs to find out the best action given some observed context (a.k.a the optimal policy reinforcement learning). Formally, we define as concretely describe the contextual
the context set and K as the number of action. Then bandit problems as follows: for each time-step t = 1,2, ... NN, where N some pre-defined time horizon (and can be given to the agent), the environment provides context s; S to the agent, then the agent should choose one action Og E 21, ... K} based on context 8t. The environment will return (stochastic) reward r(st, 0t) the agent based on the context 8, and the action o that similar to the reinforcement learning setting. And notice that, the agent can adjust the strategy at each time-step, so that this kinds of problems called "online" learning problem.
Solving the contextual bandit problems is equivalent to find some algorithms that can minimize the pseudo-regret (Bubeck et al. 2012). which is defined as:
N
N
Ry
max a:S (12,...)
E
Σ
rs,, ((x))
Σ
r(st, at)
(7)
15

---

Under review as a conference paper at ICLR 2020
Boston
Concrete
Energy
Naval
-3.04
-0.5
-3.06
4.0
-0.6
3.8
-28
3.10
3.6
-3.12
3.4
-0.8
3.2
-32
-3.16
3.0
-0.9
20
100
20
40
100
Num of samples
Num of samples
Num of samples
Num of samples
WineRed
WintWhite
Yacht
-1.04
0.4
-0.98
-0.6
-1.00
-1.06
-0.8
-102
-1.8
-1.0
-104
-12
-1.10-
11.4
1108
-166
-1.10
-1.12
-1.8
-111
-7.0
40
100
20
40
Num of samples
Num of samples
Num of samples
Figure 10: Comparison between SRLD and Langevin dynamics test log-likelihood. The results are computed based on 20 repeats. The error bar calculated based on log-likelihood of SRLD log-likelihood of Langevin dynamics in 20 repeats rule the variance of data splitting.
where ı denotes the deterministic mapping from the context set S to actions {1,2. K}} (read- ers can view π as a deterministic policy in reinforcement learning). Intuitively, this pseudo-regret measures the difference of cumulative reward between the action sequence at and the best action sequence πt (st). Thus, an algorithm that can minimize the pseudo-regret ( can also find the best и.
Posterior sampling (a.k.a. Thompson sampling: Thompson 1933) one of the classical yet success- ful algorithms that can achieve the state-of-the-art performance in practice (Chapelle & (01T). It works by first placing an user-specified prior n9. on the reward r(s, a), and each turn make decision based on the posterior distribution and update it, i.e. update the posterior distribution uáa with the observation r r( Q1-1) at time t -- where a is selected with the posterior distribution: each time, the action is selected with the following way:
Gr arg max f(st,a), f(st,a) uh.. aE (1,2,.....)
ie., greedy select the action based on the sampled reward from the posterior, thus called "Posterior Sampling". Algorithm summarizes the whole procedure of Posterior Sampling.
Algorithm Posterior sampling for contextual bandits
Input:
Prior distribution μ. time horizon N for time 1, N do
observe context 8 S. sample the reward each action f(8,, a) ~ ph.a- a E (1,2, K}. select action a, = arg max, (1.2,.... f(se. a) and get the reward r(s ax), update the posterior Mata, with r(se. 0.). end
Notice that all of the reinforcement learning problems face the exploration-exploittiii dilemma, so as the contextual bandit problem. Posterior sampling trade off the exploration and exploitation with the uncertainty provided by the posterior distribution. So if the posterior uncertainty is not estimated properly, posterior sampling will perform poorly. To see this, we over-estimate the uncertainty, explore too-much sub-optimal actions, while if under-estimate the uncertainty,
16

---

Under review as a conference paper at ICLR 2020
fail to find the optimal actions. Thus, it is good benchmark for evaluating the uncertainty provided by different inference methods.
Here, test the uncertainty provided by vanilla Langevin dynamics and Self-repulsive Langevin dynamics on two of the benchmark contextual bandit problems suggested by (Riquelme et al. 2018), called mushroom and wheel. One read (Riquelme et al. 2018) to find the detail introduction of this two contextual bandit problems. For completeness, we include it as follows:
Mushroom Mushroom bandit utilizes the data from Mushroom dataset (Schlimmer 1981 which includes different kinds of poisonous mushroom and safe mushroom with 22 attributes that can indicate whether the mushroom is poisonous or not. [Bundell al.] (2015) first introduced the mushroom bandit by designing the following reward function: eating safe mushroom will give a +5 reward, while eating poisonous mushroom will return reward +5 and -35 with equal chances. The agent also choose not to eat the mushroom, which always yield 0 reward. Same to [Riquelme et al. 2018), we use 50000 instances this problem.
Wheel To highlight the need for exploration, (Riquelme et al. 2018) designs the wheel bandit, that can control the need of exploration with some "exploration parameter" 8 E (0.1). The context set S is the unit circle ||s||2
in R², and each turn the context s; is uniformly sampled from S. K = 5 possible actions are provided: the first action yields constant reward N (N o²); the reward corresponding to other actions is determined by the provided context 8:
For E S s.t. ||s||2 I 6, all of the four other actions return suboptimal reward sampled from (2u, 0²) for μ2 < P1.
For E S s.t. || > 5, according to the quarter the context is in, one of the four actions becomes optimal. This optimal action gives reward of (( 0²) for μg pı, and another three actions still yield the suboptimal reward N (u2, o²).
Following the setting from (Riquelme et al. 2018), we set μ1 = 1.2, μ2 = 1.0, and #3 = 50.
When 5 approaches the inner circle ||I|| < 8 will dominate the unit circle and the first action becomes the optimal for most of the context. Thus, inference methods with poorly estimated uncer- tainty will continuously choose the suboptimal action 01 for all of the contexts without exploration. This phenomenon have been confirmed in (Riquelme et al.) 2018). In our experiments, as we want evaluate the quality of uncertainty provided by different methods, hard for existing inference methods as shown (Riquelme et al.) 2018), and use 50000 contexts for evaluation.
set 5 0.95, which is pretry
100
00
-05
-15
Figure 11: Visualization of the wheel bandit (5 0.95), taken from Riquelme al. 2018
Experiment Setup Following (Riquelme et al. 2018), we use a feed-forward network with two hidden layer 100 units and ReLU activation. We use the same step-size and thinning factor c/n = 100 for vanilla Langevin dynamics and SRLD, and set M = 20, 1 on both of the mushroom and wheel bandits. The update schedule is similar to (Riquelme et al. 2018), and we just change the optimization step in stochastic variational inference methods into MCMC sampler step and use a higher (2000) initial steps for burn-in like the warm-up of stochastic variational inference methods in Riquelme et al. (20T8). As this is an online posterior inference problem, only the
17

---

Under review as a conference paper at ICLR 2020
last 20 samples to give the prediction. Notice that, in the original implementation of Riquelme al. (2018), the authors only update a few steps with new observation after observing enough data, as the posterior will gradually converge to the true reward distribution and little update is needed after observing sufficient data. Similar to their implementation, after observing enough data, we only collect one new sample with the new observation each time. For SVGD, we also use 20 particles to make the comparison fair.
18

---

Under review as a conference paper at ICLR 2020
B THE DETAILED ANALYSIS OF SRLD
B. SOME ADDITIONAL NOTATION
We use ll Il+o denote the & vector norm and define the L.o norm of function f RR + R! as Il C.... Drv denote the Total Variation distance between distribution f1, P2 respectively. Also, supe, 8, ( 02). For simplicity, we may | In the appendix, also Φ (0) g(0;p), where g(0: p) is defined in the main text. For the elearance, we define TMfefn * ps: = Pk TM,e/n * pk = P and
K is R Rd R', we denote ||
||
TM. Pe == p, where M. Ph and p are defined in main text.
B.2 GEOMETRIC ERGODICITY OF SRLD
Before we start the proof of main theorems, we give the following theorem on the geometric godicity of SRLD. It is noticeable that under this assumption, the practical dynamics follows an (Me/n + order nonlinear autoregressive model when k 2 Mc/n/
0++1 2 (Ok...., On OR- Mojn) + / 2nex.
where
M
ab (Ox: O2 Me/n)) = Où + - VV(0k) + add
Σ
80 0 (0k)
M
Further, if stack the parameter by
O k
Ok Mc////
and define v (Ox)
[6TT (x)) [ we have
k+1 = V (Ox) + Ex,
where Egg - T
In this way, formulate O₂ time homogeneous Markov Chain. In the following analysis, we only analyze the second phase of SRLD given some initial stacked particles + Mc/n-1.
Theorem B.1. (Geometric Ergodicity) Under Assumption choose y, and such that
4.1 and Assumption 4.2
suppose we
max (1 1--2nas b 20m
bi, 20m (bi +11 <1,
then the Markov Chain of O₂ stationary, geometrically ergodic, i.e., fe., any O 2 Mc/n-1- have
DrV [p* (., (1), e
where some positive constant, Q(Oa) is constant related Oo. pk is the k-step Markov tran- sition kernel and II is the stationary distribution.
We defer the proof to Appendix B.5.1
B.3 MOMENT BOUND
Theorem B.2, (Moment Bound) Under assumption 4.2 suppose that we have Eo |101 < oo; and 02 a |||]] (2b1 + > O, we have
sup Eamps |10112 v sup Egnp.
|||| 2 sup Egnp, ||01|
b + +
<Bo~po po 1 2
a2 - |Klle Coo
2a a [K] (all (2b1 + 2)
And by Lemma
B. we thus have
sup Ea-pi (V) 2 V sup Egnp, |VV v sup Corp (V)22)
&
(b + 1+
<b, SbiEa-p 110|
+
02 - |||]] Co a | (2h + +
19

---

Under review as a conference paper at ICLR 2020
The proof can be found at Appendix B.5.2
B.4 TECHNICAL LEMMA
Definition B.1. (a-mixing)
For any two o-algebras A and B, the a-mixing coefficient is defined by
a)A,B)=
sup
IP -P
AEA,BEB
Let (X,, 1) be a sequence of real random variable defined on AL,,, This sequence is a-mixing if
a(n) sup (Mks (k+n) - +0, as n + 00; k2l
where Mj = (Xi,i S j) and Gj = (Xi, for 2 1. Alternatively as shown by Theorem 4.4 of Bradley (2007
Cov (f.g) f oo (MR), g E fro (Gutn)) }
a(n) = sup
[f 11g1lc..
Definition B.2. (8-mixing)
For any two o-algebras A and B, the a-mixing coefficient is defined by
B (A, = sup sup
ΣΣ
- -P())B ==1==1
where the supremum is taken over all pairs of finite partitions (A1..., A1} and { [B,,.., BJ} of Ω such that A; A and Bj E B for each i, j. Let (Xk,k k 2 1) be sequence of real random variable defined on (St, A.P)), P). This sequence is 3-mixing if
((n) sup 3 (M&, GR+n) - 0, as n co.
K21
Proposition B.1. (3-mixing implies a-mixing) For any two o-algebras A and B.
1 a (A, S 3
This proposition can be found in Equation 1.11 of Bradley 20055
Proposition B.2. A (strictly) stationary Markov Chain is geometric ergodicity if and only (()) O at least exponentially fast as n + 00.
This proposition is Theorem 3.7 of Bradley 2005)
Lemma B.1. By Assumption
4.2 we have [[V < b1 (I10:l| 1) and |0 VVVI <
(1 21ja1 + 1)²b1) |10] b 21ph1-
Lemma B.2. (Some property of RBF Kernel) For RBF kernel with bandwidth o, we have
||] 3000 < 1 and
|K(00, 01) ((0',02))
11 ///
10, 02l12
Lip
2 (1)2/7 (-)
K(O,, 01) - K 02)l
||0; 0,
Lip
Lemma B.3. (Some property of Stein Operator)
For any distribution p such that Ea-p
[V)
< oc, we have
||(()(())l < 10
Ee- VV + -2- (1//(a))
Lii
Lip
|VV(0')++ 2 (110') + (101))
[ < |Kll Earn
( + ||0'| + ||0||
< K by + Egrr
20

---

Under review as a conference paper at ICLR 2020
Lemma B.4. (Bounded Lipschitz of Stein Operator) Given 0', define Φo (0)
)(oer((0)
K (6)) + ViK(0', 8). We also denote 0g.(0)
00 . (0)
Øo a [0)]
We have
d
o (O))
s 2 2
11-11//
+ 2d
|0ll ////a 01
Lip
Lip
1
Σ
(0)]:
2d
2
1111 /0
01
+2 +2 11-111//011
VVV(0')
i=d
C00
L00
B.5 PROOF OF MAIN THEOREMS
B.5.1 PROOF OF THEOREM B.I
The proof of this theorem is by verifying the condition of Theorem 3.2 of An & Huang 1996). Suppose [0,,... O 0sc++ where = =e/n, we have
M
|/ [+)]
17
1 VV) + add
Σ
80 + (0x)
M
01 7VV())) +
70
0.133/0 (01 01t5c) -=10=,0-0:0///VVV
M
M
0 (VV) +
2 na
0 00+ 00=s0-0011//0/
M
j=1
M
+
7
10 0.11//0 2
(- (-VV(01+)c) 01t5c)
M
j=1
2017 || 1x0,00 by (1(1
$ 0|0
VV
110
M
2am
| 110,00 b1 (1+ (1+ bi |10++5oll)
+
Mo
1=1
(1)
40m
$ (1
+
|K] oo))) t1²bb1 + 271b2
σ
+ 1 21j1
+ 2017 || [K|l 100, b1 |1011| + 2007 | |Kll00000 (b1 +1)
+ r
max
|
iE[MC+1]-(1)
4am | lo) + nitb1 + 27b1
< bi (1 +
σ
+ max
1 1 + n b +
200
||||| | by
20m
||
| (br +1)
max sE[MC+1]
||01,5oll
where (1) is by Lemma B.1 Thus, given the step size n, if we choose 1, a such that
max
[ 2na, n2b4 + 2am
||
bn a
207
|| ||llo (42 +1) 1,
then our dynamics is geometric ergodic.
B.5.2 PROOF OF THEOREM B.2
Continuous- Time Mean Field Dynamics (5)
Notice that as our dynamics has two phases and the first phase can be viewed as an special case of the second phase by setting O, here only analysis the second phase. Define U. = supy
|0 and thus
Øt U? < E - v + a * p.)(O.) VO.
21

---

Under review as a conference paper at ICLR 2020
Now we bound E (0., -V(0) + aolTMM (o))):
E (0., - V(0.) + TM. * p.)(
<b1 - -2 00
oE 0 || * (0.)l
(
Sb B
0
a
K
E
0
[IV + (ell)) + +10+1l)
b1 (1100|++1) 11 2 (101) + [10/1))
R
0.
+
||
K
00
E
0:1
E
2a
2 2)
1b
02
|
co
E 0
+
K
10000
E
0:
Ee EO
b1 + o) |10'1 + b2
<h
02
K
20
U? + 0 K | R 01
Egr Bor-as
(1 (1++22)
2 |10'1 + b1
M
2a
+ 2
1
<h
|
80
U? + a |
bi
U,Ur-je + a |K| b,U,
M
==1
2a
bi + 2
<b
02
|||]
o0
U? + a |
U²2 +
||||] 1c bu? +1?
0
2a
2
S (b) +1)
02
||||||
||||]
1mx
2b +
-
U.+
0
where (1) is by have
B.3 By the assumption that x = 02
|
2a a || (2bg + 2) > 0. we
we
0 U < [bb +1 AU?] VO. 0t
By Gronwall's inequality, we have U? < U3 + bا# (If à U? O. then U. fix and this bound still
8t
holds.) Notice that in the first phase, = 0, we have X < 02 and thus this inequality also holds.
Discrete-Time Mean Field Dynamics Similarly to the analysis of the continuous- time mean
field dynamics 5 we only give proof of the second phase. Define U = supp
E 0 and thus
U2 - U S /2//EE ( (O--1. 00-1,. + abLTMMenn * pa](Ox)) ) + 21/² vo.
By a similarly analysis, we have bound
E { -VV + aofTTM * p))(
$ {b1 1) - AU?
where \ 02-
|||]
2a
|| || 11,,00 (2b) + 2) > 0. And thus we have
10,00 ø
U U2 1 < [2n] [(b] + +1) -U + 2n v 0.
It gives that
Uß < bi + 1
+ +UZ.
x
Practical Dynamics 3 The analysis of Practical Dynamics 3 is almost identical to that of the discrete-time mean field dynamics and thus is omitted here.
B.5.3 PROOF OF THEOREM 4.1 AND 5.1
Notice that the dynamics in Theorem 4.1
special case of that in Theorem 5.1 and thus we only
prove Theorem 5.T here. After some algebra, we can show that the continuity equation of dynamics
66) is
Oipt V ([- (D(0) + Q())) (V a (.)( Pt + (D(0)++(00) Vp))
22

---

Under review as a conference paper at ICLR 2020
Notice that the limiting distribution satisfies
a.d. V Q(1))))) a pxol + (000) Vpoo)
= --((((D)))Q(())
ex (xx P + + Vpx)
(1-(D(()-Q(()) VV())) + (D(0) + Q(0)) Q(0))
+av- (K + (Vp> - (Vp=>VV (O))
which implies that Poc x exp(-V(()) the stationary distribution.
B.5.4 PROOF OF THEOREM 4.2
In the later proof we ca represent the quantity
b1
+ 11 [ 100,00
EΦ-po 101
02
2a
(2b)
2)
10,000
Recall that there are two dynamics: the continuous-time mean field dynamics (5) and the discretized version discrete- mean field Dynamics (4). Notice that here couple the discrete-time mean field dynamics with the continuous-time mean field system using the same initialization. Given any T v,N, for any 0 1 S T, define [ We introduce another continuous-time interpolation dynamics:
(V))) + dB,
t [0, Mc)
0:
(V(O)) + (OM[TM, (..( dB,, 12 Mc.
P: Law (0.);
0.- On Po:
Notice that here we couples this interpolation dynamics with the same Brownian motion that of the dynamics of 0,. By the definition of è.. at any t4 key for some integrate 0 has the same distribution. Define p. = Lo,( conditioning on 0g- = 0, and p Law (
[N]. 0., and
conditioning on 0, 6g. Followed by the argument of proving Lemma 2 in Dalalyan 2017 if Mc we have
k 2
DKL
Pix
On
Ou Pie
1
vV- + aOfTTMc * (1) + VV((.) CooTTM,e p.](O.) | ds
E
4
1511
(v. - adITM, 0,, (v aO(TMM p=l(0.) ds
E
1=0
3
VV(
E
1=0
3a
t5
E B|ooTTMc (,,,( 3bTM, po](0,,) || ds
3a
k-1
1511
E BlooTM - p.1(0.))' ds
4
1=
=h Iz
Iz.
We bound /,, Io and I3 separately.
Bounding 4 and 13 By the smoothness of VV, we have
VV( vv 0 || < 03 0. - 0,
23

---

Under review as a conference paper at ICLR 2020
And by Lemma
B.3
(Lipschitz of Stein Operator), know that
1 M...M * pe)(O1) M[TM,C * pa)(O2)
< 6
1////
[V + E
2 (0)()))
10:
02112-
Lip
Iiip_
And by the Assumption 4.2
and that Ps as finite second moment, we have | p.l(0)) @bTMM p.)(O2)]
<Cca 00 02ll2.
Combine the two bounds, we have
30
2
h-1
11 13 -
d
E 0c, 0. ds.
4
=
Notice that 0, =2++ [ (0.)
aOfTM, (0) (1-1)+ Se dB.. By Itō's lemma, it implies
that
3C
2
0. 1|
4 +
13
S
E
0
ds
==
3Cc
(V + aafTMs 2=1(02))' (s (s) + 20(s - t²) ds
E
1-1
= Cc 3
Σ
E
vV(),,) + a * pe,1(0.,) )|
+ Celdkh²
j=0
By the assumption that E
0.,
is finite and 0. o E |
2 is also finite, we have
E
(V
adHTM, P2((2)
| * (0,) 2
<2R
(0.)
+ 20²E
2 + +2a2E 02
((=+b)
<467 + E
Eo Bor 110'1 + [101]
sec.
Thus we conclude that
T + I3 S Cea (cghenz + dkn)²)
Bounding Iz
MTM.C * (t,,( - @ITM, pal(O,,) 1²
M
2 plle, (1 0 (0,,))
1
M
M
plp,, (0()) - olp.-co ||?
E
M
Σ
E
Ba- 200, (0) - -E ~00, (0)
M
M
Σ
Eo,,
Ea-pay -at -0)),, -E (0,
M 1m1
iml
Eo,,
(10e) cm
v (1(ee,
Di [ārs-ot Po-cl
1=1
1=1
24

---

Under review as a conference paper at ICLR 2020
By Lemma B.4 and the Assumption 4.4 that V is at most quadratic growth and that Pr has finite second moment, we have
d
O v (1) Lip )²
0
==1
(108
v 110 (3|
-Ba,,
1.pp/
==1
4d
2
110 12/
01
+4
e 10:1///0
Bo., VV 112
BL
BL
<C(d ct)-
Plug in the above estimation, have
3a k- 1
1+
OoTMMoo ),,((O.,)) ObTM. * pal(0ss)
12
D [ōrg-ct Pa-ct] ds
5=0
S <a 5
DeL P (aa-cclDa-cll Pa-cll ds,
definition
Pinsker's
where the last inequality is due the relation that DD
DD
I BKL-
Overall Bound
Combine all the estimation, we have
t5+1
DEL
Pie
1
< <a²C(d 2
DEL
do + Cej (cght + dkm?)
M
f=0
M
a Cdd ci)
DKL
( : ( ds + Ce. (E@tnn dkn)))
M
==
1=1
Similar, if k: < Mc 1, we have
DKL.
Ptu
Ptx
0o
VV (V.)) - VV(O.) ds
E
12 k-1
E 00 0, ds
4 j=0
13113
k-1
|
dkti
Σ
E
(V,(),,)
VV
12
1=0
dkbřnn
12
4
Define
Ug
sup
DKL.
1e 1:000
00
sE
sE(te.ta+1l
25

---

Under review as a conference paper at ICLR 2020
and Us max (E( u1-
We conclude that for k 2 Me for any k s k,
"
up
S C
DKL.
- (e ds + Ce. (c2hon + dkip)
j=0
Sa
Ocs (cglin² dkm)))
mnl (1
M
k- Uj + Cia (c@kn)² + ddwr²)
sa a
Σ
For < Me, which a simpler case, we have
Ux S C ke + dkn?) < C
We bound the case when k: 2 Me.
k-
Us S a² C(d + cd)
Σ
Ug + +cd (cakn² dkn²)
1=0
If we take 77 sufficiently small, such that e
< dkjj²
we have
k- 1
UR < C(d (2))
U₁ + 2Cegdka) Σ
j=0
k-1
S o² C(d+ + 2) Σ (U, +n) j=0
Define v/ a²C(d 2n and we choose 7) small enough such that 1/ < 1/2 and 1) 1/2.
Without loss of generality, we also assume o/ 2 and thus we have
k-1
U S v Σ (U)+)))
Also we assume Us 2 r/', otherwise we conclude that UR < 1//- We thus have Up S 9 + Uy where q 2n/- Suppose that UMe M---11
CMc (nc² d) n and some algebra (which reduces to
Pascal's 's triangle) shows that
U x + 9)
We conclude that Uo 20(1 + 9)a h-1
Notice that q 2a2C(d cg),. Thus for any h: 2 Mc/n.
U S xq(1 + 9)a 4
xq(1 9) (Mn-Mo)/n
= 2q(1 + ) 202 ((d+t))(oo
S x2a² c2) e ((dt Mc),
I CMcaa (nc² d) (d + 2)) c(d+ 234+
for sufficiently small n. Combine the above two estimations, we have
Us s (nikt² + dkuz n)) CMaaa (nca + d) (d + 20 cl(+)3)(k + Cn k Mc/n
S Mc/n
Notice that now we have Uk Ie 1e(00...))sEju.n)
max sup DeL [p ll p which is function of on-
We then bound Us max sup DRL piyts Il pty]- Notice that the KL divergence has the
following variational representation:
DRL[P1 ll p2] sup Ep,, f e
26

---

Under review as a conference paper at ICLR 2020
where the f is chosen in the set that Bpt and Epg exist. And thus have
DEL|Pin+s
|| ptn) sup Eeumpo
(Egg, Pin-a
E, B,00
e!
Ptn
S Eeorpo sup
(
o0 e!
Pig-a
Pia
And thus Us & Uk. Also the inequality that
Us = max sup DKL Pinta ll ponl 2
max lE(o....k)
DKL lptn || pin)
holds naturally by definition. We complete the proof.
B.5.5 PROOF OF THEOREM 4.3
The constant h1 is defined
h1 =
2
10 2
01
||e -
1-/1//
IBL
BL
Lip
Now we start the proof. We couple the process of 0g and 02 by the same gaussian noise in every iteration and same initialization 0,- Og. For k S 1c Oz On 0 and for & 2 Mc/n we have the following inequality,
2 0x
θ OR+I
0x+1
0
22 0r 0r: -V + vVx)))
M
+2TaBBB
O2 -
0l
60
s0/////1 ((0k) 6FTM.e/n" * (o((O
M
j=1
+ t E
(V)()))
aωl
100- + VV( - aOOTMO (Oxs
M
By the log-concavity, we have
E 00 0s, -VV(OO) vV
< - LIE 02 0/
for some positive constant L. And also, as r is small, the last term on the right side of the equation small term. Thus our main target bound the second term. We decompose the second term on the left side of the equation by
E
O
O
Se S0 00 ((0x)
<|TM,o/n pol(Os)
M
M
E
0 0%
Se Sen-sern (Okk) ΦfTM,//n * px)(Ow)
M
E
0: 0R. MlTMMo// pk](0x)
o[TM.c// (0x)
+ 0: 0.. MTTMe//n (0)(0x) 6(TM.e// ph (0x)
11 1/ 13.
We bound h, Ta and 13 independently.
27

---

Under review as a conference paper at ICLR 2020
Bounding 4 By Holder's inequality,
M
1 S
02
0ull
Φ
de
30/n ((0s)
M/TMs Pk)(Os)
M
==1
M
1
Ox
0x||
E
Φ
g0 jo/n (0x)) M/TMee/n ps)(Os)
M
We bound the second term on the right side of the inequality. Define
E 0
Σド= 80 son-ss/,, ((0s) MFTM,n/n pk)(Os)
02 sup
sup |
E
M S0, se/1l(0) < MTTMM/// * pa] pa]
<B
and by the regularity assumption we know that
sup 1
E
1 M Σn 60 (00x))((0)) - & *
k
02
<0.
inf sup E |
db マ= 60, )e,,1(00 ofTM.o// px](
B
Define of M و -0 (0)
M/Tn (0) 60 1 Σド= 50 60 20/m/// and since the stein
M
operator is linear functional of the distribution, have
M Se
ESTh
(0) 0,
M
given any 0. By Theorem B.I that O is geometric ergodicity and thus is 3-mixing with exponen-
tially fast decay rate by Proposition B.2 And by Proposition [B.1] we know that O₂ also a-mixing with exponentially fast decay rate. We have the following estimation
E
50
10/4 (0k)
OfTM, (0x)
M
M
<ag sup
50 setn 100)
0 SB
M
jml
M
02
sup
E
16 160 ke/4 (0)
M
|
B
k=1
+ 02
sup
IE
16" 100, ko/n (0), --
10 1())
M
|
kテj
Co2
+
M
rnC
for some positive constant r that characterize the decay rate of a mixing. Combine this two estima- tions, we have
aagC [e-re (1
h S B 0x oall
[ 1- cra
+1
V
M
Bounding 12 By Holder's inequality, we have
1 < v E
0n /
0z
V [|o[TMM./0 * pk](Ox)
dlTM,c/n * pb)(O,))'
28

---

Under review as a conference paper at ICLR 2020
We bound the second term in the right side of the inequality.
E
MfTMM.// Pk](Ox)
M(TMMc// p=)(O>)]
M
E
plpk-jo/v|(BA) jc/k
M
E (plox- ||0(0==Sc/n((00) bc/k-)//
M
M
Σ
Eon
-. R0- (00
20 (0))
M
1=1
M
Σ
Eon
E (0) Ea- 50 (0)
M
j=1
M
Σ
(1000.(()lcc v (1 DDEL (Pk-so/n: Px-sonn
M
/=1
By Lemma B.4
we have
Σ
(1om. v (1lo)))()1)))))
d
Σ
(1
1=1
4d e
1111 //
0
+ -l 11-11-1///
(+)))
BL
BL.
Plug in the above estimation and by the relation that DBL
W1 S W2, we have
E
M,c
(0x(Ox)
MMTMM px)(0x)1 2
4d
1101
2
/00
01
+ 4
e
111///
Eo [V+())]
De. [px-ejs Pk_csl
BL
BL
M
4d
110||
/0
01
2
+ 4
-
1/11//1
Eeg
(0x)
1
M
W3 Lok- Pk-csl
BL
BL
M
And combined all the estimation and by the definition of weweeeestteiniisttnnee conclude that
1 M
1 <
4d
-
10112/00
00
+ 1111%/
2 BL
Es. (0x))
W3 [pk pk-
M
BL
1 M
4d
11010 /
| Es. [VV)) + 4 1-11e-1-111
Σ
E
Ox ej
0 0x-cs
M
BL
j=1
Bounding Is
By Holder's inequality,
Is S E
B 0x 0z
dlTM.c// * pxl(0x) ofTT * pxx(0,)
We bound the last term the right side of the inequality. By assumption and Lemma B.3 have
Φ / * (0s) - 3FFM.e/n * (xc
2
<
107//
Ea-p.
[V +
2
E
0x
0% ???
Tip
Lip
29

---

Under review as a conference paper at ICLR 2020
And combine the estimation, have
13 -
||e-0090
E Ea-pa
[VV
2
0
E 02 0z 00ll
Lip
Lip
Overall Bound
Combine all the results, we have the following bound: for h 2 Mc.
0% |
E Ox+1
Oz+1
E
0x
2
1
L1LEE
0x
+ +2na V
02
0n
2
C1
VM
M
+20002
Σ M
E
Ok-oo/n 0 02-00//
TR
0
;=1
+2nac3B
0x 02
where
re (1 e rM
C1
agC
+1
erc
2
1011 10
+ +4 +4||e--11P/00
BL sup E., (0
BL
639 /0
sup Ee-iz [|V]())) +
2 - /a)))
C3
Lip
Lip
and
C4 sup Me/n
E (V( aΦp
80 30._)e/1110) VV aølTMe/n po)(0x))
Define / E 10 0r 0a] l2
and Up = sup Ur we have [E(k)
2nac1
U2 S qU? +
Us + 17²C4,
VM
where q = (1 27((T. 0x22 - oc3)). By the assumption that a S L.((2 + (z), q 1. Now we (2naes +(1 (I=a))(c++++))
prove the bound of Up by induction. We take the hypothesis that Up <
var
(1-9)²
30

---

Under review as a conference paper at ICLR 2020
and notice that the hypothesis holds for Up 0. By the hypothesis, we have
2mar1 VM
+ (1 - 9)n
C1 + 1-9
2
2naci VM + (1 (1-9111 9)m1 (ct + 1
2nocc
+ n² c + -1a)
Uk+11 S q
(1 9)² 2
VM
(1 (1-0) q)
C1 + r )
2maass
+ (1 - 9)v7
VM
21j0xe1 VM
21j0x1 VM
- + n²
c4 + 1-9 1
+
1)
c4 + 1 1- q
(1 q)
1- q
2mac1
+ (1 - 9)n
Ca + 1-1 )
VM
(1 9)²
(1 9)²
9
21aaci
VM
+ (1 - 9)
27
VM
(
c1 +
1 1-q)
9)
+ (1 -9)1²² (1+ (a+ 1-
2mac1 vM
+ (1 - q)v
C4 + 1
<
(1 1 q)²
9
22aac
2
27ac1
VM 17 (
1 1 + +(1 10--91 -9) (cr + 1 1- q-
1 2
+ (1 q)
(at c
(1
0)²
VM
2maci
+ (1 - 9)n
c4 + I-9
2nacr vM
+ +11- - /)1l css + ( )²
VM
+ (1 9)
(1 1
(1 9 9)
2
2nac1
VM + +(1 - 9)1
C1 1-9
1
(1 9)
where the last second incquality holds by (1-9) (cr + T-1) 2 1. Thus we complete the argument of induction and we have, for any k,
2 VTT
+(1 (1-9)) 9)m (cs + 1211)
U2-
<
(1 (1-08 9) 2
2(1 - 9) 7?
(ca+ + 1
1-9
< 2
(1 2
2a
1
+ 4n7 (ca + 20(LL a022 acy))²
(L 002 ac2)² M
And it implies that Wilpk, Po] su S Un S
(L (L-ac2-aca) M + 417 (cs + 217(L - ac2 - -ocz))2
20
B.6
PROOF OF TECHNICAL LEMMAS
B.6.1 PROOF OF LEMMA B.
For the first part:
[0)
< (0) -V0)
<b [1]
31

---

Under review as a conference paper at ICLR 2020
For the second part:
10 VV((
-VV
= 1101 + 2n (0,-V((0) + + (0V
S |10l² + 2n (aa |012 + b1 + 1 |011²)
(1 2qa + +/bbb) 1101² +,1²b 2npr
B.6.2 PROOF OF LEMMA B.2
It is obvious that | ll s0 S 1.
[K(0', 01) K(0',02)
100'-00
01|l2///// 112/0 -c-10--0a 111//1/
< ()) 110, -22122
Lip
And
Vo K K(' 01) 1.i-oo Vo-K(0',0
- 1100'0011 -111'/0 /a (0' - -0,) -,))- - -100''0221' (0' 0z)
< - 0 (/())) 110 - 091-22112
B.6.3 PROOF OF LEMMA B.3
For any distribution p such that (O) < oo,
(0,) - (12)
Ea-p {- [K(0, 01) - ((0000)) -
(1)2/a
Eap [V] [| - -2112
Lip
2 - (1) (a())))))))))))))())))))))))((
||0: 02112-
Lip
For proving the second result, we notice that
[()p Egp (((O,,)VVV(O)) V,K(0'00) 0)]
|VV 2 ([10') ([ll])
< || K Eg'op
(b+b) | |8
< || K 3o b1 + EBrrp
B.6.4 PROOF OF LEMMA B.4
Given any 0',
32

---

Under review as a conference paper at ICLR 2020
d
00 ,s (0)
i=1
(01) 60 (02)
sup
|00 2 l
00
d
0 (01) 200,,, (02)
sup
0, #0
|0 02l12
c 1|0' 0, 11"
-- ||0' 100'-0/l|//0
a V(0))
2
00;
sup
01 /02
110 -2212
|10r -0.10
11. 11 / (0 - 01)
2 e
e' -021 (02.s 07)
+2
sup 01 v00
||0| 02 2
For the first term on the right side of the inequality.
(c-110' 10:10'-0:11-- 1//1//
-ollo'-Bal1%/0 V(0') 2
80
sup
01102
|0 02112
(e~10'' 0|0 0, 11 // -e 10'-0al1/10
i
2
Σ
V(0")
sup
09,
0,
10, 0211
(0)
11-11
//
Lip
To bound the second term, by the symmetry of each coordinates, we have
2e -0000 11.1 //0 (01, 03) - 10000000 (0, -00) 2
sup
#z
00 0
2
0113 1013//0000 01
d
Lip
This finishes the first part of the lemma.
d
j0 (0))
Σ
-1 0' 01 -////
2
0:
2
0
V(0')
V
2
0
O0
d
d
2 e 000 01
(0, 8)
+ Σ 2
100
00 v
Amd
6
imd
c..
2
2 000
01 /0 (0: 0¹)
12
+2 1.11/%
[VV(
i=d
Cx0
<2d
2 e
|1011
/
+2 e
111////
C.
Coc
33


Reviews:
[{"review_id": "rJeF90yCKS", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper proposed another variant of Langevin dynamics, called “Stein self-repulsive dynamics,” which simultaneously decreases the auto-correlation of Langevin dynamics and eliminates the need for running parallel chains in SVGD. They combined Langevin dynamics with Stein variational gradient descent and theoretically justified that the proposed method successfully converges to the stationary distribution with only a single chain, unlike SVGD. The proposed method decreases the auto-correlation of Langevin dynamics, so the proposed method increases the sample efficiency.\n\nThe paper is well-written. The idea of the proposed method is natural, which is incorporating the functionality of SVGD to reduce the auto-correlation of Langevin dynamics. The idea is intuitive and justified by their theoretical analysis. The authors also well- placed their work in the literature, as described in Section 3. The intuitive explanation of the proposed method is given in Section 3.\n\nI have one technical question as follows. If the authors reply appropriately, I will raise the score to accept.\n\nIn Theorem 4.3, the result holds for any k and M. The authors claim that if we take a limit of M -> ∞ with fixed k, the practical dynamics converges to the discrete-time mean-field limit, in Section 4. However, to state the result of Theorem 4.3, k should be bigger than M c_\\eta from the dentition of \\tilde{\\rho}_k^M, as shown under the equation (4). How do we take a limit of M -> ∞ ? Does k also go ∞?\n\nMinor comments:\n- The definition of g should depend on only \\theta_k^I and \\hat{\\delta}_k^M, not \\theta_k^k.\n- The equation (1) should hold for any \\theta’, not \\theta.\n- The equation (1) should contain \\rho, not p.", "extra_scores": {"experience_assessment": "I have published one or two papers in this area.", "review_assessment:_thoroughness_in_paper_reading": "I read the paper at least twice and used my best judgement in assessing the paper.", "review_assessment:_checking_correctness_of_experiments": "I assessed the sensibility of the experiments.", "review_assessment:_checking_correctness_of_derivations_and_theory": "I assessed the sensibility of the derivations and theory."}, "extra_text": {}}, {"review_id": "H1xeyoxG9S", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "I have the rebuttal of the authors, the paper improved indeed and some point on role of M is better clarified now although it is still a bit convoluted. The paper would be stronger if the analysis shows any theoretical advantage to the presented method. I think the author put a good effort in addressing some of my concerns and I m raising my score to 6. \n\n \n####\nSummary of the paper:\n\nThe paper proposes stein self repulsive dynamics for sampling from an unnormalized distribution. The method starts by using Langevin dynamic for up to time $Mc$ and then uses those pasts samples to guide the trajectories of the langevin sampling to explore new areas of the densities  using the stein witness function between the current particles  and the past samples( similar to Stein Variational gradient descent). \n\nThe paper analyses the mixing properties under standard assumptions of the potential of the Boltzmann distribution and the kernel used in the Stein discrepancy, and shows convergence to the boltzman distribution as the number of particles goes to infinity and the step size goes to zero. \n\nAuthors validate their methods and show that it indeed explores on a synthetic example new areas wrt to pure langevin dynamics. Applications in sampling from the posterior of bayesian neural networks and contextual bandits compare the performance of the proposed method to langevin dynamic and pure stein descent favorably. \n\nClarity/presentation :\n\nThe paper is well written and the intuition are well presented . \n\nThe notation $\\hat{delta}_{M} $ is not great in denoting direct measures, please using another symbol. \n\nMy main concerns with the paper are the following:\n\n- the definition of $\\bar{\\delta}_{M}$ averages only $M$ particles choosen for $M$ time stamps. Something is off here for the continuous approximation to work at each time stamp you need $N$ particles and then you have a past horizon $M$. As $M\\to \\infty$ this does not matter, but I think in your implementation you are considering at each time step $N$ particles and you average on a horizon of size $M$. is this correct? Please clarify?\n\n- The theorem show only asymptotic behavior and don't quantify the intuition behind the paper , that the \"coverage of the samples\" is higher.\n Can you for instance bound the wasserstein  distance between the pure langevin and Stein Repulsive dynamic , and between SVGD and your method, as it was done in \"The promises and pitfalls of Stochastic Gradient Langevin Dynamics\". Basically you can find a coupling between trajectories of your methods and the langevin dynamics and bound the wasserstein distance between the two methods. This will be insightful to see if one would mix faster with respect to the other one.  \n\nWhile the appendices of the paper are lengthy I don't think they explain the most selling point of the method, since they are asymptotic and there is a confusion between the time horizon for the past , and the number of particles.  if $M$ is  $\\infty$ in this defintion , we are just running langevin dynamic, and not using Stein witness function. I see an important issue in this definition of the time horizon, I think this time horizon should be finite, and that from time steps less we sample N particles , and we let this $N$ go to infinity.\n", "extra_scores": {"experience_assessment": "I have published one or two papers in this area.", "review_assessment:_checking_correctness_of_experiments": "I carefully checked the experiments.", "review_assessment:_thoroughness_in_paper_reading": "I read the paper thoroughly.", "review_assessment:_checking_correctness_of_derivations_and_theory": "I assessed the sensibility of the derivations and theory."}, "extra_text": {}}, {"review_id": "BkgV-wG5FH", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The papers described how to use the repulsive term used with standard SVGD within MCMC/SGLD. Briefly, the paper proposes to use a (damped version of) the SVG repulsive term between the current position of a SGLD trajectory and the empirical distribution defined by the trajectory.\n\nThe approach is interesting and natural. Unfortunately, I do not think that the experiments are convincing. \n\n(1) Since the authors are advertising the Bayesian framework, the choices of metrics such as test RMSE or test LL are not adapted\n(2) in the UCI dataset examples, it is indeed extremely difficult to explain the enhanced performances? Is it because of multimodality? Better exploration of a mode? Comparison to a single NN? Comparison with an ensemble of NN? \n(3) I would have been much more convinced by a set of well-chosen and controlled experiments. The 2-dimensional examples are far too low-dimensional to be convincing. Higher-dimensional Gaussian? Higher-dimensional mixtures? Non-linear tractable problems in higher dimensions? Influence of the tuning parameters (RBF parameter? alpha? step-size in SGLD)? Computational issues? Subsampling-effect? etc...\n\nThe proposed method is interesting and has a lot of potential. I would like to suggest the authors to spend more times on careful and controlled numerical experiments (Bayesian NN are not very good for this purpose) -- with convincing numerics (which would give more reasons to delve into the proofs) the method can be very promising.\n", "extra_scores": {"experience_assessment": "I have published in this field for several years.", "review_assessment:_checking_correctness_of_derivations_and_theory": "I did not assess the derivations or theory.", "review_assessment:_checking_correctness_of_experiments": "I carefully checked the experiments.", "review_assessment:_thoroughness_in_paper_reading": "I read the paper at least twice and used my best judgement in assessing the paper."}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposes a new sampling mechanism which uses a self-repulsive term to increase the diversity of the samples.

The reviewers had concerns, most of which were addressed in the rebuttal. Unfortunately, none of the reviewers genuinely championed the paper. Since there were a lot of good submissions this year, we had to make decisions on the borderline papers and this lack of full support means that this submission will be rejected.

I highly encourage you to keep updating the manuscript and to rebusmit it to a later conference.
Decision: Reject