Path Auxiliary Proposal for MCMC in Discrete Space
Authors: Haoran Sun
Hanjun Dai
Wei Xia
Arun Ramamurthy
Venue: ICLR 2022 Spotlight (2022)

Abstract: Energy-based Model (EBM) offers a powerful approach for modeling discrete structure, but both inference and learning of EBM are hard as it involves sampling from discrete distributions. Recent work shows Markov Chain Monte Carlo (MCMC) with the informed proposal is a powerful tool for such sampling. However, an informed proposal only allows local updates as it requires evaluating all energy changes in the neighborhood.
In this work, we present a path auxiliary algorithm that uses a composition of local moves to efficiently explore large neighborhoods. We also give a fast version of our algorithm that only queries the evaluation of energy function twice for each proposal via linearization of the energy function. Empirically, we show that our path auxiliary algorithms considerably outperform other generic samplers on various discrete models for sampling, inference, and learning. Our method can also be used to train deep EBMs for high-dimensional discrete data.

Published as a conference paper at ICLR 2022
PATH AUXILIARY PROPOSAL FOR MCMC IN DIS- CRETE SPACE
Haoran Sun
Hanjun Dai Google Brain hada @goog le com
Wei Xia"
Georgia Institute of Technology hsun3496ga edu
Amazon
xxi aww.iaaaaazzor
Arun Ramamurthy Siemens arun. www.mmmttaaryssimesssn.. com
ABSTRACT
Energy-based Models (EBMs) offer a powerful approach for modeling discrete structure, but both inference and learning of EBM are hard as it involves sam- pling from discrete distributions. Recent work shows Markov Chain Monte Carlo (MCMC) with the informed proposal is powerful tool for such sampling. How- ever, an informed proposal only allows local updates as it requires evaluating all energy changes in the neighborhood. In this work, present path auxiliary al- gorithm that uses a composition borhoods. We also give a fast version of our algorithm that only queries the eval- uation of energy function twice for each proposal via linearization of the energy function. Empirically,
local moves to efficiently explore large neigh-
show that our path auxiliary algorithms considerably outperform other generic samplers on various discrete models for sampling, in- ference, and learning. Our method can also be used train deep EBMs for high dimensional discrete data.
INTRODUCTION
Many real-world problems involve discrete structured data modeling, such as syntax trees for nat- ural language processing Tai et al.) (2015). graphical model for molecules( et al. (017). etc. powerful approach for modelling the distribution over structured data Energy Based Mod- els(LeCun et al.) 2006) (EBMs). EBMs define the distribution with an unnormalized energy func- tion, which allows great flexibility to fit the target distribution. However, this ficxibility also results in the difficulties in inference and learning, as they require sampling from the EBM (Andrieu al.
2003
Hinton 2002), where the partition function is intractable many cases.
Markov Chain Monte Carlo (MCMC) algorithms are one of the most widely used methodologies to sample from intractable distributions (Robert & Casella 2013). The efficiency for MCMC depends drastically the proposal distribution. For example, in continuous space, Metropolis-adjusted Langevin algorithm (MALA) exploits the gradient of the target single step walk and biases the proposal distribution toward high probability region (Rossky al., 1978 Roberts & Rosenthal] 1998 Welling & Tch. 2011); Hamiltonian Monte Carlo (HMC) employs multi-step walk and ex- plores the distribution more efficiently(Neal [2004] Girolami & Calderhead 2011 Hoffman al.) 2014). These methods substantially improve the performance of the MCMC algorithm in theory and practice. However, their proposal distributions are derived discretization of continuous diffusion process and it is still not clear how to appropriately extend such methods into discrete space(Zanella 2020).
proposed general framework called pointwise informed proposal (PIP) that shows promising results on directly sampling from discrete distributions. PIP utilizes the energy change in the neighborhood of current state to propose a new state. Following this work, Grathwohl et al. [022] propose a more efficient sampler that uses Taylor series to estimate the energy change
Recently, Zanella 2020)
WWor done during the time Siemens

---

Published as a conference paper at ICLR 2022
in the neighborhood. However, both methods only focus on proposing from small neighborhood, for example, a H-Hamming ball. This is due to the computational expense of evaluating the energy change or approximation error of Taylor series for large neighborhood. As result, the samples in Markov chain will have strong correlation or even be trapped at local optimum, which deteriorates the sampling efficiency.
In order to efficiently explore a larger neighborhood, we propose a path auxiliary sampler, which is an auxiliary sampler (Liang et al. (011) that uses auxiliary path to propose new states. In construc- tion of the path, we employ a local proposal distribution to make a sequence of small movements to reach a new state in long distance. A typical challenge for such multi-step proposal that the accept ratio could decrease very fast when the number of steps increases. In this work, we show that the accept ratio of our path auxiliary proposal independent of the path and only determined by the property of the current state and the proposed state. As result, our algorithm is able to maintain a high accept ratio with multiple local steps. We can provably show that the locally balancing func- tions are asymptotically optimal for our path auxiliary sampler. We also introduce scalable version of the algorithm that uses linearization and can be applied to smooth target distributions.
We empirically evaluate our methods in inference, sampling, and learning various discrete EBMs. We demonstrate that our methods significantly improve the sampling efficiency on parity model, weighted permutation model, Ising model, Restricted Boltzmann Machine, and factorial Hidden Markov Model. Our method can also learn competitive deep EBMs on discrete image data. The code can be found at htttps:///gt wtwpth/w.cccoooarraasnnaaaaaaaaa
ler git
BACKGROUND
Energy Based Model: Let X' be a finite state space. An EBM defines an energy function f
X' x R with target distribution (x) = e (w) /2, where Z = =xeex f(=) (LeCun et al. 2006 Wainwright & Jordan, 2008 Du & Mordatch 2019). The unnormalized energy function provides great flexibility to characterize complex distribution. However. this can also make the partition function Z intractable to calculate and exacerbate the difficulties of learning and inference.
Petropolis-Hastinges MCMC is a commonly used methodology to sample from an intractable distribution. Metropolis- Hastinng (MH) is one of the most commonly used framework for 1970). Given the current state 2, proposal distribu- {(y)Q( tion Q(x,,) gives state y. then MH algorithm accept y with probability of min x())((z,)) to guarantee the Markov chain is reveessiibe.. The efficiency of MH algorithm highly depends on the selection of proposal distribution Q.
MCMCoMetropolis et al.] 1953] Hastings
Peskun Ordering: Peskun ordering provides method to compare the efficiency for two MCMC algorithms (Peskun 1973 Tiemey 1998). Let P1, P2 be r-reversible Markov transition kernels on X' such that Pi (a y) > c 2 (z. y) for all # y for a fixed e > O, then we say P1 c-times more efficient than P2 the following holds: var, (h, Pi) < { var, (h, P2) v= (h), Yh : X > R; 2) Gip(Pi) 2 e G2p(P2). A smaller asymptotic variance var, (h, P) means a better estimation for the expectation of and larger spectral gap Gap(P) P)) means a faster convergence of the Markov chain.
Pointwise Informed Proposals and Locally Balancing Function: A pointwise informed pro- posal (PIP) is MH algorithm in discrete space (Zanella 2020). PIP uses proposal distribution Qg(x, y) = o( m(y) x())) [Ix ()) where the [(x,8) v) = (vEN(()) is the membership indicator
w.r. symmetric neighborhood Ni(-), and the normalizer (z(r) Σ EEe 9 T(()) [ 2). The
x(E)
weight function g + + R determines the efficiency of PIP Zanella (200 show that the family
of locally balancing functions g {g R4 - RR.g(e g(t tg Vt O is asymptotically opti- mal for PIP w.r.t. Peskun ordering. Empirically, Zanella 2020) shows g(t) vi and g()) z have the best performance.
Gibbs with Gradient GWG (Grathwohl et al. 2021) is a scalable version of PIP, where the target distribution (--) is approximated via Taylor's series. Despite being powerful, PIP and GWG requires the calculate the weight and sample from the neighborhood. Hence, only a small neighborhood, usually 1-Hamming ball, is used in existing methods.
2

---

Published as a conference paper at ICLR 2022
3 PATH AUXILIARY PROPOSAL
Problems occur in point-wise informed proposal (PIP) when merely considering small neighbor- hood, especially for distributions with many local optima. For example, consider a parity distribution where the state space X {0, 1}" and the energy function f(z) U if the number of Is in is odd, otherwise f(x) = O. When a 1-Hamming ball neighborhood is used, by symmetry, PIP will have a uniform probability to propose a new state from neighborhood. Then, the expected time to leave low energy state is O(e"), which could be very inefficient when U is large. When r--aamiing ball larger neighborhood is used, PIP can efficiently escape a low energy state. How- the neighborhood will contain Onn") states, which could be computationally prohibitive when n or is large.
To address such problems with PIP, propose path auxiliary sampler. Instead of directly sam- pling from large neighborhood, path auxiliary sampler sequentially samples new state from a local proposal distribution Qn- The computation cost each step is still manageable Qo still samples locally. Then the composition small movements forms path that can lead to a new state that is distant from the current state. The complexity of the sampling grows linearly instead of exponen- tially w.r.t. r. In this section, we first present the framework of our algorithm. Then, we discuss the choice of the weight function for path auxiliary sampler. And finally, we introduce fast version of path auxiliary sampler.
3.1 PROPOSAL VIA AUXILIARY PATH
We first define an auxiliary path. Given the neighborhood function N, the set of auxiliary paths on X' is defined
[( ((,,
L) Eii
Li N((a-))), I,.., L}
(I)
To obtain the auxiliary path, we employ a PIP Qo(:.) .) to make local movements. We also sample path length from prior a(-) to assure our chain aperiodic. The path auxiliary sampler is defined as follows
1. Sample path length L from prior (L).
Denote the current state I; = oo, sample o; ~ Q0(01-1,.) for I....,L.
3. Accept T++1 = L in probability A(x, o, L) = min {1, )(का) IT (o(m) { else F+41 = Dt- Theorem 1.
त( l1 2ola 1.71)
The path auxiliary proposal transition rule described above satisfies the detailed bal- ance and induces a reversible Markov chain with π as its invariant distribution.
3.2 BENEFIT FROM AUXILIARY PATH
We first analyze the performance path auxiliary sampler on parity distribution mentioned above. Assume a uniform prior (<1) a(2) 1.2 If we sample L = 2, then will transit back to a new state with same energy as the current state. As a result, we can always escape from state in O(1) steps in expectation.
Having path auxiliary sampler to sample from L-Hamming ball at each step can be more efficient than performing L steps of MH sampling in 1-Hamming ball. To mathematically justify it, compare their accept ratio. For path (a, L):
min 1, T(OL) IL Color, Co(ov,oL-)) M (TO) IhE Qu(ar-1,0))
II min {1, ra(i)Q(la
1=1
We can notice that the accept ratio for path auxiliary sampler is the product of the probability for single-step method before the minimum operator. As result, auxiliary path can have high accept ratio as long as the product large. On the contrary, when performing L steps MH sampling, the transition probably be blocked some steps in the path with low accept ratios.
3

---

Published as a conference paper at ICLR 2022
3.3 BALANCED PROPOSAL
Although the product in equation alloow the Markov chain to escape from a local optimum when some steps in the path have high accept ratio, it also exponentially decrease the accept ratio w.r.t the path length L when every step has low accept ratio. Hence, is important to select a good local proposal Qo. In this section, we show that the locally balancing function Zanella| (2020) is asymptotically optimal for path auxiliary sampler. We first define a sub-class of weight function named as ideal function:
Definition 1. F is the set of ideal function. A function f E F if and only if following conditions
holds: 1) f R 1R+ 2) [(1) = 3) [(() monotonic increasing:4) [())(()) S 1,
The next theorem shows that locally balancing function is asymptotically better than ideal function F in Peskun Ordering.
Theorem 2. Consider the state space in Cartesian products X x"= where each X, is a finite space with M elements, and the neighborhood is defined as I-Hamming ball. Let d. be the maximum degree in conditional independence graph. If 1) lim du 0; 2) the larget distribution satisfies m(y) SC oc, Vy E Ex(x; 3) the path length prior a is bounded by U. Then for any f F. we ((())
have function g E -= (9()) gt) tg(()) that is asymptotically more efficient than which implies g is asymptotically optimal in F.
The idea to prove this theorem to use f(t) = V f(t)f(t). By definition, we have:
f(t) f(t)= = - N
(3)
which means f E G. Then, the theorem is proved by showing the x-reversible Markov transition kernel Py(x, 3) 2 Py(x, y) asymptotically holds. The proof is given in Sec A.5
Although Definition hhas some constraints, it already contains almost all natural choices of weight functions, such as (f) i++ f(t) min f), ft) max(1, f). and f(t) = to, 2 O. Hence, the asymptotic optimality for locally balanced function G in ideal function F. strongly suggest that locally balanced function is a good choice for path auxiliary sampler. More discussion for ideal function can be found in Sec[A.3]
2t
Property 1. When using a locally balancing function g G as the weight function, the accept ratio for path auxiliary sampler can be writlen as:
A(x, V, min { Zg Zz(x)
(4)
7,(y)
The path independent form in Property shows that if the energy function is smooth enough, the accept ratio in equation 4 will be high when y is close to x. By selecting g(t) = vi E G, we name auxiliary sampler (Liang al. 2011) ) as path auxiliary sampler (PAS) and give the algorithm in Algorithm
3.4 SCALABLE PATH AUXILIARY ALGORITHM
Although path auxiliary proposal reduces the complexity from exponential to linear, it still requires ((nL) evaluations of the probability <xx) in total. In cases where the computational bottleneck is evaluating Txx). this cost will be expensive. Fortunately, most distributions of interest have differ- entiable energy functions, such as deep EBMs, and we can use the linearization as approximation (Grttwohl et al. 2021). Given the current state de, the linearization is
(2) f(ao) + (Vf(oo), -ou)
(5)
estimate No. (v)) (I: Q (x, = o( Ro(y)) ) ) /2000 (I), with the normalizer Zy,cro (a) = EveN(w) 9( Fo(w) (r) After sam-
e (V)()))))))) and define our proposal distribution fo(z)
pling o, we use the linearization Jor at or to compute the accept ratio w.r.t. the auxiliary path
A(r, L) min T(v) IT Qg.n. (01,
(6)
((r) I (0L-1,
4

---

Published as a conference paper at ICLR 2022
Algorithm 1: Path Auxiliary Sampler (PAS) and the fast version (PAFS)
Input: Target Distribution Initial state zo, Path Length Prior a, Weight Function g(t) Vt Output: Sample sequence (I1, T2)..)
repeat
Sample path length L ~ a(L)
Denote do 2t, sample or Qo(or-11)) or ~ Quin (01-1.) for L.,,.
Accept 7++1 0r. in probability A(x, a, L) = min(1, 2z(0l)
[o(we)
or A(z, L) = =0=mm] min {I, x(y) IT Cu,y (0. else T+41 Tt.
xxx) [T= 10111))
until finish sampling;
By using approximation, we only need to evaluate (v)), π (x), Vi()), (zx) once in every pro- posal, which can significantly reduce the computational cost. We name this algorithm as path auxil- iary fast sampler (PAFS) (sce Algorithm|| D
The PAS framework also allows different approximations. For example, GWG Grathwohl et al, 2021 can be seen as special case where all indices are sampled from Q (oo:) rather than Qg.oo ().. Besides, one way to further encourage long range movements is sampling sites without replacement, that's say, the auxiliary path does not modify a site more than once. However, such a sampler requires a decent choice of the path length L. When L is too large, the acceptance rate drops exponentially fast. Hence, in this work, we focus on PAFS. Similar to GWG (Grathwohl) et al. 2021), the decrease of the proposal quality from PAFS can be bounded.
Theorem 3. Assume the energy function (x) is differentiable, Vf(x) is K-Lipschitz. Consider we use weight function g(l) Vt, path length prior a bounded by U, and I-Hamming ball as neighborhood. Denote and P as the transition kernel by path auxiliary sampler and path auxiliary faster sampler, respectively, have
P(x,y) y) 2
P(x, y)
Using the Peskun ordering, we know PAFS is at least e e -K200 200++11 2 times as efficient as PAS. When K small, theorem 3 justify PAFS good approximation for PAS. When is large, tighter bound needs more assumptions for the target distribution, e.g. conditional independence. See the proof and discussion in A.7
RELATED WORKS
Informed proposal for Metropolis-Hasting (MH) algorithm has been extensively studied in the
continuous space(Robert & Casella 2013). The most famous algorithms are Metropolis-adjustél Langevein algorithm (MALA) (Roberts & Rosenthal 1998) and Hamiltonian Monte Carlo (HMC). MALA, HMC. and their variants (Girolami & Calderhead 2011 Hoffman et al. 2014 Welling & Tch 2011 Titsias & Papaspiliopoulos 2018) exploit the gradient of the target distribution to bias the proposal distribution towards high probability regions. Although gradient-basee baseed methods hav- ing brought substantial improvements continuous space, is still unclear how to extend them to discrete space.
number of methods try to map the discrete space continuous space using relaxation, apply informed methods in continuous space, and then map the new state back into discrete space(Zhang al. 2012, Pakman & Paninski] 2013) Nishimura al.) 2017: Han Liu 2018: 8 Jaini et al.) 202 VAE flow. Such methods work in some scenarios, but a key challenge is the embedding of discrete space into continuous space can destroy the inherent discrete structure, resulting in highly multi-modal and irregular target distribution in continuous spaces.
via Gaussian Integral Trick, uniform dequantization,
Another group of methods directly work discrete space. Dai et al. 22200) introduces the path as latent variable in the variational distribution for initializing PCD, but still relies on slow Gibbs sam- pling for improvement: Titsias & Yau (201)) augment the discrete space with auxiliary variable and
5

---

Published as a conference paper at ICLR 2022
perform Gibbs sampling in the augmented space based on informed proposal. Zanella (2020) shows that a family of locally balancing function is asymptotically optimal for informed proposal. Follow- ing Zanella (2020). Power & Goldman (2019) extends the framework to Markov jumping process, and (Sansone (022)) parameterize locally balanced function to tune it via mutual information to se- lect good weight function from the locally balanced class. When the target distribution is smooth enough, Grathwohl et al. (2021) employs Taylor's series to approximate the target distribution and further improve the sampling efficiency. Though these informed proposal algorithms successfully show orders of magnitude improvements when using 1-Hamming ball neighborhood, they are not able to explore large neighborhood. tially address this problem via sampling multiple dimensions to modify in one step. However, such procedure can easily lead to backtracks thereby reducing its efficiency.
extension of GWG (Grathwohl al.) (202T) can par-
EXPERIMENTS
5.1 SAMPLERS UNDER CONSIDERATION
In this section, we empirically evaluate the sampling efficiency of path auxiliary sampler (PAS) and path auxiliary fast sampler (PAFS). We choose the path length prior a as uniform distribution on {1....,XX 2X- }} as suggested in Hoffman et al. (2021) and denote the corresponding samplers as PAS-X and PAFS-X.
We compare our methods with five types of baselines; random walk sampler (RW), Gibbs sampler (Gibbs), Hamming ball sampler (HB), locally balanced sampler (LB), Gibbs with gradient sampler (GWG). RW is informed proposal with (t) = 1 that uniformly propose new state from neighbor- hood. Gibbs partitions the dimension of a state r into two groups Iu and x-u- then updates the state from conditional distribution p(zrul>-.))- We denote as Gibbs-X, where X refers to the dimension
of Tu- HB is a two-stage Gibbs sampler on the extended state space (x, x²) (Titsias & Yau 2017
We use HB-10-1, where 10 is the block size, and 1 is the hamming ball size. LB (Zanella 2020)
implemented version of LB, which draw the index to flip via first order Taylor's 's series of the target distribution. Although the original paper focuses on flipping one site per step, it is possible for GWG to draw
stggle-tepp PIP with gt) VĖ. GWG (Grathwohl et al. 2021
scalable
multiple indices to flip in every step. We denote the algorithm as GWG-X, where X indicates that
the number of indices to modify per step is uniformly sampled from (L.... 2X 1}. The difference between PAFS and GWG is that GWG is more likely to sample the same index repeatedly thereby reducing the efficiency especially after mixing, e.g figure
5.2 INFERENCE ON ENERGY BASED MODEL
Parity Model: We first demonstrate the benefit of path auxiliary sampler versus single-step sampler in inference task that estimates the mean f of a parity distribution. A parity distribution has state space X' {0, 1 3 and energy function
f(x) = (Ei_)) 2x) mod 2)U
iee. a state has energy U if it has an odd number of 1s, otherwise O. The neighborhood defined as 1- Hamming ball. We run the simulation with p 100 and U =(1, 3, 5}. and report the estimation
100
LB
2000 usoa 1oox 2000 32.00 30804 Energy fumiian Endlester
J8.3 1500 X814 2500 3000 Dvagy Tencuen Dedoodono
5K0 13000 3000 2500 J00N 4X000 Energe Fundoon Evenattrr
Figure 1: Estimation Error of Distribution Mean on Parity Model
error En : || e where Ar, = E'= X, is the sample mean of the Markov chain step n. For each setting and method, we run chains for 20,000 steps, and we plot the mean and standard deviation of the estimation error in figure T] We can observe that the efficiency of single step sampler

---

Published as a conference paper at ICLR 2022
100, Sgma 33
p. Sgmm
-10
LB-1
RW-3
-220
RW--3
LB-1
-50
-00
LB-11
RW-33
-30
IM5-3
PAS-3
PA5-5
-50
.12
2000
421
2000
15904
2823
5000
5002
1090
2000
4000
8000
18050
Figure 2: Energy Trace on Weighted Permutation Model
therey
0443
3000
MCN: S1831
MO3C Samp
Figure 3: PAS and PAFS in different lengths on 200 x 200 Ising Model
decreases exponentially when U increases. On the contrary, the path auxiliary sampler always escape from the local optima in parity distribution and estimate the distribution mean efficiently.
Weighted Permutation Model: We consider optimization task weighted permutation model. The state space X' = Sp is a symmetric group, i.e. for any P E X, p is a permutation of 1, 2, ..., p. The energy function is defined as:
(p)== Ei- ((()(()
Following Zanella (2020), the weight {wi Fijj= are i..d. sampled from (aussian(0,00)) The neighborhood determined by one exchange of the permutation. We use p 100, o 10, and run 5 simulations for each configurations. From figure 2 we can see that IB sampler is trapped local optimum, especially when a is large and the distribution is sharp. On the contrary, PAS can always efficiently find good solutions via auxiliary path.
5.3 SAMPLING ON ENERGY BASED MODEL.
Lattice Ising Model: Consider the state space X' = 131, where (Vp, Ep) is the x p square lattice graph. For each x E X, the energy function is defined as:
O, ~ (L,EE,, 213,
(10)
0g R are bias terms representing the property of an, and is global interaction term. We first compare PAS and PAFS in different path lengths. For each length, we run 100 chains with 100,000 steps. We report the burn-in stage as well as effective sample size (ESS) in terms of both MCMC steps and energy function evaluation times in 3 We can see that, T) PAS has slight improvement in terms of energy function evaluations when we increase the path length; 2) when compared to PAS, PAFS has very similar proposal quality and significantly better efficiency. We also compare our sampler with other competitors. For each sampler, we run 5 Markov chains with 1,000,000 steps and report the ESS to compare the proposal quality in figure 4 sampler leads in both quality and efficiency. More results are given in See B.
We can see our
factorial Hidden Markov Model: FHMM is statistical model that use latent variables in X {0, 1 N to characterize time series data y E R.W. Denote p(r) for hidden variables, and p(y(z) for likelihood:
NN II
x II p( II p(em,clomn-)) k=
p(u]>)
Gaussian (Yn; wr, + b, 0²) ==1
(11)
Given data y. we use MCMC to sample x from the posterior p(rly) and compare the mixing for different samplers. We choose parameters N 1000, K 10, P(x,,1 = 1) 0.05,
'Computed follows Tensorflow Probability

---

Published as a conference paper at ICLR 2022
Evolustions
of Reconstruction
107
Ga0-1
540.
107
107 50P 35
03
-05
10F
-1.
-2.
Energy Function Evaluations
Figure 4: ESS for Different Samplers on 200 x 200 Ising model
Figure 5: Burln FHMM
P(Xn.k = In.k-1) 0.85, w E Gaussian(O, Ix). ). E Gaussian 1), and o² 0.25. We each chain 5 times and report the mean and std for the energy (negative log joint density) E(ē) log ((v))( and the reconstruction error |y ly-wat will2 w.r.t. the number of evaluations of energy function. From figure we
that both GWG and PAFS gain significant improvements
in mixing by using larger neighborhood. More results are given in Sec B.2
Restricted Boltzmann Machine: RBM is bipartite latent-variable model, whose energy function is defined as:
f(z) = =og(1 + e==x+o) + 6²
(12)
where (W, b. c) parameters and {0,1)P We follow Grathwohl et al. 20021 train RBM with 500 hidden units on the MNIST dataset using contrastive divergence(Hinton 2002
We generate samples via various MCMC samplers the trained RBM. Besides reporting ESS. we also estimate the maximum mean discrepancy (MMD) (Gretton et al. 2012) w.r.t. a set of "ground truth" samples generated by structure known Block- Bibck sampler. Figure 7 shows PAFS has large ESS and can efficiently match the "ground truth" samples. We also demonstrate an ablation study on different the approximations of the auxiliary path. Specifically, we evaluate GWG, PAFS 100 chains and report their ESS and average hop distance. We can notice GWG suffers from a large path length L, while our PAFS obtains robust
with different lengths. For each sampler,
1.0
GAN
GAG
600
GMG
6X4
300
2. 302 2. Average Length
3.3 74. 25.3
125
124 Aucrago langth
oge Path Length
Figure 6: Sampling RBM with Different Path Length
improvements for all L by employing soft no-replacement sampling. This result indicates that our PAFS provides an efficient framework to sample from large neighborhood. More results and discussions, including no-replacement sampler, given in Sec B.3
RBM logIMMDD)
RBM ESS
17"
-159
GWC
-135
Gw2-48
OMSO
Gw2-
M956
-22
19
NUSSS
-1.
NUS-5
Tnery/ Funoson Tasliniiis
Simpling Stopa
Figure 7: Sampling on RBM trained on MNIST
Figure 8: Learning Ising
8

---

Published as a conference paper at ICLR 2022
VAE
VAE
EBM
EBM
EBM
Dataset
RBM
DBN
(MLP)
(Conv)
(GWG)
(Gibbs)
(PAFS)
Statie MNIST
-86.05
-82.41
-80.01
-117.17
-79.58
-86.39
-85.67
Dynamic MNIST
-82.42
-80.40
-80.51
-121.19
-79.59
Omniglot
-103,52
-97.65
-94.72
-142.06
-90.75
-100.47
-100.78
Caltech Silhouettes
-112.08
-106.35
-96.20
-163.50
-84.56
Table 1: Evaluation of different discrete models on the held-out test set.
5.4 LEARNING ON ENERGY BASED MODEL
Learning an EBM is challenge task. Consider the target distribution is x and our energy function fe is parameterized by The gradient for the likelihood of TT(x)) x e -f(()) is:
V log (z) = E,, fofo(()) I oo
(13)
The first expectation can be estimated using the data from true distribution. The second expectation requires samples from the current model, which is typically obtained via MCMC. Hence, the success of training an EBM relies on efficient MCMC algorithms to get an accurate estimation of the second expectation.
Ising model We first learn an Ising model following the setting in Grathwohl et al. (2021 where the energy function f(x) 0x Ja. Given set of samples {2²): the task to learn an EBM via recovering the connectivity matrix J. We generate 25 x 25 2D cyclic lattice with 0 = 0.25 and sample training data with a long-run Gibbs chain. We train the models to maximize the likelihood of samples using persistent contrastive divergence
(Tellmaa & Hinton 2009) and report the root mean squared error (RMSE) between the inferred connectivity matrix J and the true matrix J. For fair comparison, we also include GWG-eq, which runs GWG-1 with more steps such that has the same time cost as PAFS-5. Figure shows larger neighborhood can effectively help learn the Ising model. GWG-5 having similar performance as PAFS-5 as, using PCD, the learning efficiency is mainly determined by the mixing of the Markov chain.
Deep EBM We evaluate path auxiliary sampler by learning a deep EBM. The experiment follows the setting in (Grathwohl al. 2021). We train deep EBMs paramterized by Residual Networks(He
iet al. 2016) on small binary image datasets using PCD(Tieleman & Hinton 2009) with replay buffer Du & Mordatch 2019). We compare our methods with Variational Autoencoders Kingma &
Welling (2013), GWG (Grathwohl et al. (021). RBM and Deep Belief Network
Hetwork(Hintoo (Hinton] 2009). We estimate the likelihoods using Annealed Importance Sampling(Nes) Neall (001) and report the results in table I The results for VAE is taken from Tomezak & Welling (2018). RBM and DBN are taken
from Burda et al. (2015), GWG taken from Grathwohl et al. We can see that our approach
improves the log- lieliiiissoods for deep EBMs on all datasets.
6 DISCUSSION AND CONCLUSION
The problem for the selection of the path length remains open. In this work, we propose a soft no-replacement method PAFS. PAFS is not sensitive the path length larger than the optimal length, hence has robust proposal quality for different choice of path length. We also considered no-backtrack samplers that flip sites without replacement. Intuitively, such auxiliary path encour- ages to propose states in larger distance. However, for no- no-backtrack samplers, a too large path length harmful for the proposal quality and significant reduce the efficiency of the sampler. As a result, finding a principle way to decide the path length is very important, In continuous space, the optimal step size can be characterized via acceptance rate (Gelman et al. 1997 Roberts Rosenthal] 1998] 2001 Beskos al. 2003).. PAS is the gradient based sampler in discrete space, we believe it is possible to derive its optimal path length in a similar manner as continuous case. We will investigate it in our future work.
In summary, informed proposal has shown good results for inference, sampling, and learning EBMs in discrete spaces. Our path auxiliary sampler provides an approach allowing informed proposal to efficiently explore large neighborhoods. We believe there is considerable room for future works to improve the sampling methods in discrete space.

---

Published as a conference paper at ICLR 2022
Acknowledgement: This research was supported in part by the Defense Advanced Research Projects Agency (DARPA) under Contract FA8750-20-C-0542 (Systemic Generative Engineering). The views, opinions, and/or findings expressed are those of the author(s) and should not be in- terpreted as representing the official views or policies of the Department of Defense or the U.S. Government.
REFERENCES
Christophe Andrieu, Nando De Freitas, Arnaud Doucet, and Michael I Jordan. An introduction to meme for machine learning. Machine learning, 30(1):5-43, 2003.
Alexandros Beskos, Natesh Pillai, Gareth Roberts, Jesus-Maria Sanz-Serna, and Andrew Stuart. Optimal tuning of the hybrid monte carlo algorithm. Bernoulli, 19(5A):551-1154, 2013.
Yuri Burda, Roger Grosse, and Ruslan Salakhutdinov. Accurate and conservative estimates of mrf log-likelihood using reverse annealing. In Artificial Intelligence and Statistics, pp. 102-110. PMLR, 2015.
Hanjun Dai, Rishabh Singh, Bo Dai, Charles Sutton, and Dale Schuurmans. Learning discrete energy-based models via auxiliary-variabll local exploration. arXiv preprint arXiv:2011.05363, 2020.
Yilun Du and Igor Mordatch. Implicit generation and generalization in energy-based -ased models. arXiv preprint arXiv:I903.0866 2019.
Andrew Gelman, Walter R Gilks, and Gareth O Roberts. Weak convergence and optimal scaling of random walk metropolis algorithms. The annals of applied probability, 7(1):110-12) 120, 1997.
Justin Gilmer, Samuel Schoenholz, Patrick F Riley, Oriol Vinyals, and George E Dahl. Neural message passing for quantum chemistry. In International conference on machine learning, pp.
PMLR, 2017.
Mark Girolami and Ben Calderhead. Riemann manifold langevin and hamiltonian monte carlo methods. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73(2):
123-214, 2011.
Will Grathwohl, Kevin Swersky, Milad Hashemi, David Duvenaud, and Chris J Maddison. Oops i took a gradient: Scalable sampling for discrete distributions. arXiv preprint arXiv:2102.04509, 2021.
Arthur Gretton, Karsten M Borgwardt, Malte Rasch, Bernhard Schölkopf, and Alexander Smola. A kernel two-sample test. The Journal of Machine Learning Research, 13(1):723-777 2012.
Jun Han and Qiang Liu. Stein variational gradient descent without gradient. In International Con- ference on Machine Learning, pp. 1900-1908. PMLR, 2018.
W Keith Hastings. Monte carlo sampling methods using markov chains and their applications. 1970.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recog- nition. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp- 770-778, 2016.
Geoffrey E Hinton. Training products of experts by minimizing contrastive divergence. Neural computation, 14(8) 2002.
Geoffrey E Hinton. Deep belief networks. Scholarpedia, 4(5):5947, 2009.
Matthew Hoffmman, Alexey Radul, and Pavel Sountsov. An adaptive-meme scheme for setting tra- jectory lengths in hamiltonian monte carlo. In International Conference Artificial Intelligence and Statistics, pp. 3907-3915. PMLR, 2021.
Matthew D Hoffman, Andrew Gelman, et al. The no-u-turn sampler: adaptively setting path lengths in hamiltonian monte carlo. J. Mach. Learn. Res., 15(1):1533-1623 2014.
10

---

Published as a conference paper at ICLR 2022
Priyank Jaini, Didrik Nielsen, and Max Welling. Sampling in combinatorial spaces with survae flow augmented meme. In International Conference on Artificial Intelligence and Statistics, pp. 1349-3357. PMLR, 2021.
Diederik Kingma and Max Welling. arXiv:1312.61 14, 2013.
Auto-encoding variational bayes.
arXiv preprint
Yann LeCun, Sumit Chopra, Raia Hadsell, M Ranzato, and F Huang. A tutorial on energy-based learning. Predicting structured data, 1(0), 2006.
Faming Liang, Chuanhai Liu, and Raymond Carroll. Advanced Markov chain Monte Carlo methods: learning from past samples, volume 714. John Wiley & Sons, 2011.
Nicholas Metropolis, Arianna W Rosenbluth, Marshall N Rosenbluth, Augusta Teller, and Edward Teller. Equation of state calculations by fast computing machines. The journal of chemical physics, 21(6):1087-1092, 1953.
Radford M Neal. Annealed importance sampling. Statistics and computing, 11(2):125-139) 2001.
Radford M Neal. Improving asymptotic variance of meme estimators: Non-reversible chains better. arXiv preprint math/0407281, 2004.
Akihiko Nishimura, David Dunson, and Jianfeng Lu. Discontinuous hamiltonian monte carlo for sampling discrete parameters. arXiv preprint
700500 853, 2017.
Ari Pakman and Liam Paninski. Auxiliary- variable exact hamiltonian monte carlo samplers for binary distributions. arXiv preprint (xXiv:/3/1.21600 2013.
Peter H Peskun. Optimum monte-carlo sampling using markov chains. Biometrika, 6063):677-612, 1973.
Samuel Power and Jacob Vorstrup Goldman. Accelerated sampling on discrete spaces with non- reversible markov processes. arXiv preprint arXiv: 2019.
Christian Robert and George Casella. Monte Carlo statistical methods. Springer Science & Business Media, 2013.
Gareth O Roberts and Jeffrey S Rosenthal. Optimal scaling of discrete approximations to langevin diffusions. Journal of the Royal Statistical Society: Series (Statistical Methodology). 60(1): 255-268, 1998.
Gareth O Roberts and Jeffrey S Rosenthal. Optimal scaling for various metropollishastings algo- rithms. Statistical science, 16(4):351-367, 2001.
Peter Rossky, JD Doll, and HL Friedman. Brownian dynamics as smart monte carlo simulation. The Journal of Chemical Physics, (9(10):422844633 1978.
Emanuele Sansone. Lsb: Local self-balancing meme in discrete spaces. arXiv preprint
arXiv:2109.03867, 2021.
Kai Sheng Tai, Richard Socher, and Christopher D Manning. Improved semantic representations from tree-structured long short-term memory networks. arXiv preprint 103.001503.00075 2015.
Tijmen Tieleman and Geoffrey Hinton. Using fast weights improve persistent contrastive di-
vergence. In Proceedings of the 26th annual international conference on machine learning, pp- 1033-1040, 2009.
Luke Tiemey. note on metropolis- kernels for general state spaces. Annals of applied probability, pp. 1-9, 1998.
Michalis K Titsias and Omiros Papaspiliopoulos. Auxiliary gradient-based sampling algorithms. Journal of the Royal Statistical Society: Series B (Statistical Methodology). 80(4):749-767, 2018.
Michalis K. Titsias and Christopher Yau. The hamming ball sampler. Journal of the American Statistical Association, 1598-111158-11111 2017.

---

Published as a conference paper at ICLR 2022
Jakub Tomezak and Max Welling. Vae with a vampprior. In International Conference Artificial Intelligence and Statistics, pp. 2214-1223, PMLR, 2018.
Martin Wainwright and Michael Irwin Jordan. Graphical models, exponential families, and vari- ational inference. Now Publishers Inc, 2008.
Max Welling and Yee W Teh. Bayesian learning via stochastic gradient langevin dynamics. In Proceedings of the 28th international conference machine learning (ICML-II), pp. 681-688. Citeseer, 2011.
Joe Whittaker. Graphical models in applied multivariate statistics. Wiley Publishing, 2009.
Giacomo Zanella. Informed proposals for local memc in discrete spaces. Journal of the American Statistical Association, 115(530):552-86 2020.
Yichuan Zhang, Zoubin Ghahramani, Amos J Storkey, and Charles Sutton. Continuous relaxations for discrete hamiltonian monte carlo. Advances in Neural Information Processing Systems, 25: 3194-3202, 2012.
12

---

Published as a conference paper at ICLR 2022
A PROOFS
A.1 PROOF FOR THEOREM 1
Our path auxiliary sampler is proposal distribution augmented auxiliary sampler. Before diving into the proof, we first introduce some background knowledge for auxiliary sampler, based on the results from page 86 in Liang al. (2011). Given proposal distribution T(ylæ:), augment it by an auxiliary variable u, such that ( = (ulx))2)222 u)du. Using this proposal distribution, we define a marginal sampler and auxiliary sampler.
1. For a marginal sampler, the accept rate is
Amar (z. y) = min(1, 7(y) (ull))) u)du,
(())
min
(14)
(2x) f (2) u)du
xxx)[(y|=)
and for y o x, its transition kernel is:
Anar (z. v) (() (x,
(15)
For an auxiliary sampler, given x, we first sample auxiliary a Ti (u(z)) then we sample the new state y ~ 2, u). The corresponding accept rate is:
x (v)TT(uyy (uly)
aux (I, y, u) min(1.
(16)
π ( (u
Since the proposal and the accept rate depend on auxiliary u, we need to integrate the auxiliary variable to obtain the transition from x to y::
K (x, y) = /. Ti (ula))aalyl,, A), u)du
(17)
Though a marginal sampler is Peskun better than an auxiliary sampler (Titsias Papaspiliopoulos 2018), marginal sampler is intractable in most of scenarios, as it requires integrate auxiliary variables u. in both proposal and accept rate calculation. Hence, our sampler is implemented in an auxiliary fashion.
Going back our theorem, the auxiliary variable u is the auxiliary path (a, L). Condition on L, we have:
Ti (ula,, L) = Ti = II Co(ar-1,a))
(18)
T2(yl x u, L) = I(x, y are the two ends of path u(or sayo))
(19)
Given path length L, with little abuse of notation, the transition kernel for our path auxiliary sampler K K(2, L) is:
L)
A(a, all)
Te(ul2.
(v) II Qolar, O1)HHo x}
Σ
II 1=1
(I(or-
H1r - y} min 1, ग (x: II 07)0 = =
()) 00=8
Tiaak, LyTa or L)A(s, n. Liddo
K(x, L) = () L)A(,, o(
(20)
When 0L * y, the accept rate above is not well-defined. To avoid making extra definition, absorb T2, the indicator, as constraint of the domain in the integration. Then, we have:
Ti
Mz.y.a/))
L (
T ( l' Qo(o).0-11) n (x:) & (0(01-1,01)
Σ
(o(o--10))) min 1=1
(21)
Og PD=3,01=y
K(,. J, Ta L)Ta(v o, L)A(s, v, L)oo With these knowledge prepared, begin proof theorem 1.
13

---

Published as a conference paper at ICLR 2022
Proof. Denote (xx, y)
the probability that transit to different state y, then have:
T ((x))((n
(22)
{1, ( I1 Qo(or,
= (< Σ
(L)
E
min {1. ((x) I C0 (o
(23)
V=1
oo=2,oL=y
min {
*(v) II Qo -
a(ll)
Σ
<(x)
II (oloo--1,0))
1, ( <(x) II Qo Q0(0--1,0))
(24)
1
(OLLEEE 00=2,01.3
1=1
- Ca(l)
Σ
min T(x)
II Qolor, 0olor,on11)
(25)
(a, N)
01
[ min
x x Co(or-1.01)
(L)
E
T(V)
(26)
T TT)) Qo(or 01-1)
(E) do 00=2,01=1 my
(1 [[Qooo)00-))
(xx Qo(or, (T) I Qo(or-1,0))
Ty)) Σ a())
E
min
(27)
() 00
TV) (v,
(28)
The key idea for the proof is that (o, L) is symmetric w.r.t its two ends, when o is a path from y. is also a path from y a. Hence, able exchange the orientation of the path in equation 27
A.2 PROOF FOR LEMMA 1
Lemma 1. is a finite state space with distribution я and a neighborhood function N such that every state has equal number of states. F is a class of weight functions, such that for any f F, f) RR+ > Rt, 2) f(1) = 3) f(t) is monotonically increasing, 4) f(t)(()) & 1, S 1. Vf E F. define f(1) - ((()(()), have
functiooo
min Zy(x)
max Zz(y) yEX
(29)
Proof. Since X is finite, always find state 3T, 22 having highest and lowest probability, respectively. Then, we have:
21(x1) = Σ || ((3
S T(2)
Σ
(z(22)
(30)
T(a1)
T
=EN(x1)
=EV(wz)
A.3 DISCUSSION FOR IDEAL FUNCTION
We define ideal function as stepping stone to show the ad advantage of locally balanced function. We name "ideal" as the properties are ideal assumption for weight function in PIP. Conditions 1) f: R4 -> R4 is natural requirement for weight function. Condition 2) f(1) easily realized by substituting f(t by f())/(()) as PIP only depends on the ratio of weights. condition 3) f(t) is monotonically increasing, and condition 4) f(t)s(t) < 1, Vt < 1, though, are technical requirements, is reasonable for weight functions. Monotonic increasing in 3) indicates the proposal distribution match the target distribution. where point has higher probability in target distribution should also has a higher probability in proposal distribution. Condition 4) is weaker than the following condition: ʃ(1) 1, vt 1, hence condition 4) is easier satisfy.
1 can be
14

---

Published as a conference paper at ICLR 2022
Currently definition of ideal function class already contains most of commonly used weight function, 2t
such as f(t) = 1+1 f(t) = max(1,t),f(t = min(1,t), and f(t) = t², where a 2 0. Also, we
can following properties
Property 2. The logarithm of ideal function is closed in convex combination. That's to say, Di=1 X L, we have
Vfi..... fn E F. and V\i,.... "a......., An S.I X 2 0 and
F(t) = e log 11(1)
(31)
Proof. Let f; and X; be defined above. Obviously, F is positive and F(1) hence satisfies condition 1) and 2). For condition 3), we have:
F(t)F()) e d f:(1) Eis Xe log fi(:) Cinm Ar log
(32)
e Ae(log (i
(33)
Eis ds (1oo
(34)
Hence, Vt S 1, we have:
F(t)S(I)t Ein 2.0 =1
(35)
which implies F(I)= F.
Property 3. Ideal function is superset for normalized monotonically increasing locally balanced
function G₁ == [g E (1) = 1, is monotonically increasing).
Proof. Ag G1, condition 1), and 3) are obviously satisfied. For condition 4), since g(t) is locally balanced, we have:
a(0at -1 (
(36)
implies E F.
Considering the ideal function cover such a large number of functions, we believe the asymptotic optimality over ideal function strongly suggests locally balanced function good choice for path auxiliary sampler.
A.4 PROOF FOR LEMMA 2
The (undirected) conditional independence graph corresponds to Markov random fields:
Definition 2. The conditional independence graph X is the undirected graph G (K, F), where K = [l....,k} and (i,j) & E if and only if X; 1 Xj1Xxetis)
More details can be found in p.60 Whittaker (2099).
Lemma 2. Consider the state space in Cartesian products X X;= where each X, a finite space with M elements, and the neighborhood is defined as I-Hamming ball. Let d, be the maximum degree in conditional independence graph. Define dim) SuPpEN(()) Z,(x)
2,1g) If1) lim
0:
2) the target distribution satisfies
x(y) 7(z)
< Yy E (xx). then we have:
du
(n
1 1 (n) 9 < 1
as n n Vg FF
(37)
(g(2))
To simplify the notation, we assume M = 2 in this proof. It is straightforward to extend the proof to any finite M.
Proof. On one side, by definition, we can easily see that og 2 1. On the other side, given I, index represent that select Ax)) by flipping index i for We de- note (i) = g(T(())) where obtained by flipping index j for æ, and we denote ge(i) = (( where is obtained by flipping index j for y. We also denote B(i)
15

---

Published as a conference paper at ICLR 2022
the Markov boundary for i, which means given B()), i is independent with remaining nodes. Then, we can write
j=11 9:(1)
d()
sup sup =EX (E[M]
(38)
j=1 9z (z)
JEB(i) 9:(0) + Σ j$B(1) 9:(1)
sup sup x ic[M]
(39)
Σ ((B(s) Io(i) + + jEB(s 9z(5)
(EB()) 9i(i)
Ejgn() 9:(i)
< sup sup =EA' ie[M]
2 jE Bi) Ii( + Σ LjEn()) Io(j) + EigB(r) Ix)
(40)
Ejen(,
< sup sup 1+ Σjen( 0n (n)() (n) E)nn() 9n
(41)
xEX' (c[M
(n.
S sup sup + nt( zE ic(M]
(42)
1 dng(C)
(43)
ng(())
Bii) is not empty, the first term in equation 40 equals to 1. If B(i is empty, the first term in equation 40 does not exist, hence it can still be bounded by
Remark: For fixed ideal function g F, Lemma 2 shows cl converges at rate + O
A.5 PROOF FOR THEOREM 2
We prove the theorem by using Lemma T and Lemma
Proof. We first show that, for all E F, the transition probability:
1 CgCg
(a( 2
(44)
where cy is defined in Lemma 2 We temporarily ignore the superscription and will add it back at the end of the proof. Consider path (a, L), denote 4, -(,))* Then the probability that do transits to y = σg is:
x (v) IT' Qy 20(0),0-12)
Po(z,8,o|D) = II Qg(0111,07) min
x (x) I1 (g(01-1,0)) 1.
(45)
[=1
g(((-1,) (g(
g(6(,1-1)
min
4--111
(46)
Zg(
By definition of Cg. we have:
min IT o(ti-1,1), II (1-119(tu-11)
P, cf
(47)
max (z)
Po(x, min [12,,0411 I l (19
(48)
mingex(0) (2)
Use the property of g. we have
L I
)
min
(th .) II (o
S
1 II
((())1)))
(49)
1=1
1=1
16

---

Published as a conference paper at ICLR 2022
Combining equation 47
equation
48
equation
49 and Lemma 1 we have:
I 1(1-11))
Po( 2
(50)
of min (z)
min [Il [11-1 I 4 (-119((1111))
(51)
L
max (2)
9
Po(a,y,o(L)
(52)
c 9
The inequality holds for arbitrary auxiliary path (o, L) and the path length & U, hence we prove the first step. Then Vg FF, we the estimation of cg
(n) in Lemma 2 we have
Pg(2,8) 2 P (x)
(53)
where
c5 = 1 OUU dng(C) ng(
(54)
This indicates I is asymptotically better than and proves the theorem.
In the proof, equation 54 shows the convergence rate depends g. this factor prevents a uniform convergence rate for the ideal function class F. We can obtain a uniform rate when we constraint in smaller function class, for example, if we restrict to FM == {f F f(C)/f())
Remark:
M). the convergence rate 0 d.
A.6 PROOF FOR PROPERTY
Proof.
( # Ooa 0-11
A(x, U, L.) = min
(55)
TTx) I ((01-1,11)
((σ))
IT o( a(o--1)
x(ar)
)/2g(o))
min
1,
(56)
((o-1))
I I x(o,)
1/2z(o1-1)
(
g( (11)))))
L (g(0)-1)
min
1,
(57)
(
g( ((()))
x(0r)
Zg(ar)
7g(x) [
min
(58)
Zy())
A.7 PROOF FOR THEOREMS 3
The idea to proof theorem 3 is similar to Grathwohl et al. (2021). which uses the Lipschitz condition to bound the estimation enror. The different is, in PAFS, the linearization at each state to propose the next state dt+1, hence to need to accumulate the error at each step to obtain our final result.
17

---

Published as a conference paper at ICLR 2022
Proof. For any path (o, L), first bound the estimation enror fu(z) fu(or-1), for E M((r-)))
Since f is M -smooth, we have:
f(2) f(or-1)
(59)
K I|= 0111
=(V((o--)),=- - -0)11) + 2
(60)
K
=(VS(oo), = - o1-1) + (Vf(or-1) - Vf(oo), - - 01-1) +
212 01
(61)
2
=(Vf(o0).2=00 + K(I )
(62)
Similarly, we also have:
k(I+ 1 1
f(z) [(o-
(63)
With this estimation, we have:
Q0(01-1,01)
(64)
Σ ExEN(O e e
(65)
ExEN(on [1(1)))
= Q(oi-1,0iee-- -Ki
(66)
Similarly, we can also obtain
@1(0) 0()1
(67)
Now, consider the transition kernel
II Qo(01-1,01) min {
x I (1(0)
P(x.v)= = Ea(L)
Σ
(68)
x(x) I (0(01-1,0))
(a () do mo=m,0==y
Ea(L)
Σ
min
(ly) Qo(or-1,1))
II (Q(0)))01))
(69)
(τx)
Go=E,OLEy
Kl (1) ((x)
2 Σ
a(L)
Σ min II Q()-
Q(a 1)- ()) -K(L-))
(70)
2 Σ
(L)
Σ
+1
min
1(00- 01), I |=1
[ly) x(x)
II Q(o))))))
b
())
00=0,00=
(71)
K 11
P(a,y)
(72)
L
Remark: Though the theorem proved above is very loose, provides framework prove ap-
proximation bound for path auxiliary sampler via single step approximation in equation 66 More
sophisticated bounds can be obtained by improving the single step estimation in equation 66 via the
property of the conditional independence. Specifically, local modifying of the state usually does not influence most of the indices, which for example can be characterized by the Markov boundary. Then, with high probability, the auxiliary path can avoid manipulating correlated indices and hence result in the approximation error in equation 66 be controlled by a constant, rather than e -Ki For example, in Ising model, an auxiliary path has high probability to avoid manipulating two adjacent nodes, and in this case, the linearization accurate and we loose nothing in approximation. In deep EBMs, the conditional independence is usually not available, but can still expect low correla- tion between most of dimensions. Such analysis relies on more accurate assumption of conditional independence and we will leave it to future work.
18

---

Published as a conference paper at ICLR 2022
65L11
+4%
14%
2445
10-3
4017
400
300
020
+800
Figure 9: PAS and PAFS in different lengths on 50 x 50 Ising Model
baalastion
-00.9
2440
-0000
-300
240-0
Figure 10: PAS and PAFS in different lengths on 100 x 100 Ising Model
B DETAILS FOR EXPERIMENTS
B.1 SAMPLING ON ISING
Following Zanella 2020), we set the interaction term À 1, and we 04 = J + Z, if in the center of the square lattice, and a; = - + 2, otherwise, where = 2 and 2; ~ 3n, 3). We generate the lattice Ising model in four sizes: p = 50, 100, 150, 200. We first run PAS and PATS with path lengh L = 1, 2. 3. 5, 10 on Ising model with size p 50, 100, 150, 200. For each path length and each model size, we run 100 chains 100,000 steps and compute the ESS using the last 50,000 steps. The results for p- 200 is given in figure 3 We give the results for p 50, 100, 150 in the following. From figure 9 figure 10 figure IT and figure 3 we can see that the difference between PAS and PAFS in terms of MCMC steps is decreasing when the model size is increasing.
The reason is that when the model becomes larger, most of the nodes are conditionally independent.
When the auxiliary path does not involve two adjacent nodes, the approximation has zero error. This observation provides a possibility to further improve the approximation bound in Theorem 3 We will study this in our future work.
For comparing with other samplers, we report the results for each problem size p = 50, 100, 150, 200 in figure 12 figure 13 figure 14] figure Our path auxiliary samplers substantially outperform other competitors.
B.2 SAMPLING ON FHMM
We compare PAS and PAFS with different path length on FHMM. We choose parameters K = 10, P( 1) 0.05, P(Xn,k 1n, 0.85, w E Gaussian (0, Ix), b E Gaussian 1), and 0² = 0.25. We simulate PAS, PAFS with path length L = 1, 3, 5, 10 on FHMM with N = 500, 1000, 1500, 2000. For each configuration, we run 100 chains and report the burn in period, as well as the ESS in terms of MCMC steps and energy function evaluations.
Damiin in terma of INOHC
ES5 Cerr NONCC HAE1
Zuren in terma of Onargy Cvalsutton
In terma ol Energy Emiletors
26-2
36-3
06/
400
245.1
64%
1250 09
MCM:
706
MOA Semsinn
350 40g 110 1996 J0hn wede
MOPC Sanstes
Figure 11: PAS and PAFS in different lengths on 150 x 150 Ising Model
19

---

Published as conference paper at ICLR 2022
ESS in terms of Energy Function Evaluations
ESS in terms of MCMC steps
Giobs1
Gibte-5
Gibbs5
Rik-
RM6-1
10a
10a
10a
10
4
6
Figure 12: Sampling on 50 Ising model
ESS
Energy Function Evaluations
ESS
MCMC steps
Gitt-1
Gibbc1
10a
GitteSS
RMk
RW-10
30
(8-)
LB-1
102
GWIS
GWYS 10
G1WS 10
PAS 10--33
PAS 13
FASF-10
PASP-1D
101
10a
Figure 13: Sampling on 100 x 100 Ising model
We can observe that: 1) terms MCMC steps, PAFS has very similar burn in time in energy compared PAS, while the the gap in ESS is larger than the gap in Ising. The reason FHMM not as close to linear model as Ising and PAFS has larger estimation error. 2) in terms energy function evaluations, PAFS still leads the performance other models. PAS only obtains faster mixing, and no improvements in ESS. The results show that whether PAS can help depends on the property of the target distribution.
B.3 SAMPLING ON RBM
We compare GWG, PAFS, and NB with different path lengths on RBM trained on MNIST dataset. GWG (Grathwohl al. 2021) samples all indices manipulate based the linearization at starting state the path dn. and we classify it as sampling with replacement. Our PAFS samples the indices based on the linearization at current state ar and we classify as sampling with soft replacement. Specifically, consider a binary case, when we flip index at state on, then in the remaining of the path, the energy change for index will be reversed. Hence, if flip an index and reduce the system's energy, we will have small probability flip it back. A third type sampler is No Backtrack (NB) sampler. also path auxiliary sampler and it rules out the indices have been selected assure
ESS in terms of Energy Function Evaluations
ESS in terms of MCMC steps
10a
10a
4
6
Figure 14: Sampling on 150 x 150 Ising model
20

---

Published as a conference paper at ICLR 2022
BLATE CANTA MOLLL Caps
Recomttration
Barrers bao wolkarat
-20
122
MCPC
bregy wwp.aaalaa
NCN: karpors
Figure 15: PAS and PAFS in different path lengths on FHMM N=500
3015
0145
0415
55
385.2
14000
04112
100
Figure 16: PAS and PAFS in different path lengths FHMM N=1000
it is sampling with hard no replacement. Specifically, NB first calculate the probability for indices based on (an) as GWG. Then, once selecting index σT, NB manipulates the probability for selecting index i in the remaining path. For each sampler with path length L = l..., 20, we run 100 chains and report the MMD (Gretton et al.) 2012) w.r.t set of "ground truth" samples generated by structure known Block-Gibie sampler. We can see that, when increasing the path length, PAFS has steady fast mixing, GWG decreases the efficiency in mixing, and NB fails in mixing. For this reason, we only compare the statistics for GWG and PAFS in our main text in figure
B.4 RUNNING TIME LEARNING ISING
We follow the experiments in Grathwohl et al. (021). GWG runs faster 1.4x than our PAFS-5 as the energy function computation is cheap this experiment. To obtain fair comparison, we add EWG--quivalentt-enge
(GWG-eq), which we allow GWG-1 to run 1.4x more sampling steps such that its running time is the same as PAFS. When our PAFS runs 5, 10, 25, 50 steps, GWG-eq runs 6, 13 34, 69 steps. Hence, when we draw figure g the point for GWG-eq at step 5, 10, 25, 50 using the RMSE obtained setp 6, 13, 34, 69.
B.5 DETAILS IN LEARNING DEEP EBMS
DataSet
Static MNIST
Omniglot
Caltech
PAFS-3
0.685
0.885
0.893
PAFS-5
0.746
0.980
0.962
PAFS-7
0.775
1.030
1.010
GWG
0.877
1.130
1.140
Table 2: Run Time for One Batch for different discrete models.
Recuntiation
Fm1
$017
1001
$15/3
201-2
$1h-1
410
001-
0700
-259
-150
RDFC 2321
SaYpars
Crergy Ciclu7a7
CATS
Figure 17: PAS and PAFS different path lengths on FHMM N=1500
21

---

Published as a conference paper at ICLR 2022
BaTeR
Racosst actor SiTa
2001
605
-20
100
40.9
2331
Sampurs
Crery Enllurters
POMC 67
Figure 18: PAS and PAFS different path lengths on FHMM N=2000
14h terns MLMC caps Loba3
40h tames MLMC Copl Lin3
B. oi teets REIIE deea L-3502
54.
2440g
60%-4
04%
34% Up
500 terre Kragg CwMakeeo Uodeo
522 turre Keagg Colaadee Le
555 Cerms Crergy www. 2533)
5 Samm creragy Eracstirn 204
2001
$003
6013
2407
Figure 19: ESS for PAS and PAFS in different path lengths FHMM
We use the same setting as Grathwohl et al. (02I), including the batch size, number of iterations, the PCD hyper parameters, etc. In principle, both PAFS and GWG requires two energy function evaluations per MCMC step, despite that PAFS explores much larger neighborhood space via path auxiliary. As the energy function evaluation and gradient calculation usually dominates the compu- tation, we expect these two should have similar runtime.
In practice, we report the average running time per batch of 100 chains in table 2 We can see for deep EBMs, our PAFS with different path lengths run slightly faster than GWG. Given that both of the methods implemented using pytorch the minor speed up might be due to the implementation issue. So overall we think it is fair to say the two methods run equally fast in most cases. For this reason, use the same number of steps when Grathwohl et al.] (202T), we the same number of steps (40) to train the binary deep EBMs. For our method, we also tune the expected path length in {3. 5, 7}, and report the result with the best validation likelihood.
train the deep EBMs. As is reported in
0.4
224
2.4 -14
-30
Figure 20: GWG Burn In
Figure 21: PAFS Burn In
Figure 22: NB Burn In
22


Reviews:
[{"review_id": "n67BxIpdh5Y", "reviewer": "Reviewer_CSos", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "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": "Strengths:\n\nThis paper presents a simple improvement over LB and GWG which alleviates the largest issue with these two samplers -- their inability to escape local energy wells. For LB/GWG to escape such a well we would need a long series of *accepted* transitions through low-likelihood regions, meaning we need the proposal to choose these unlikely transitions *and* the A/R step to accept them. Conversely, for PAS-X we only need to choose these unlikely transitions and the accept will then be likely. This removes a large factor that keeps these transitions from happening. \n\nOn top of this, I enjoy the fact that the implementation of this algorithm is nearly trivial on top of LB or GWG. One must simply add the path-length proposal and apply the LB/GWG proposal X times. Methods like this are likely to find quick impact since they can be easily incorporated into existing systems. This is compared to continuous relaxation-based discrete sampling methods which add many additional hyper-parameters, must be carefully tuned, and still tend not to perform as well as the methods discussed in this work.\n\nI find the results presented here compelling. On a wide array of distributions, PAS (and PAFS) appear to outperform LB (and GWG) when evaluated in distinct ways. Further, there *appears* (see comments in Weaknesses section) to also be a notable improvement when using PAFS to train EBMs. \n\nThe work is not overloaded with theory (a good thing!). The authors provide what is necessary for understanding and correctness of the algorithm while leaving the work accessible to readers who are not experts in MCMC. Theorems 1 and 2 were easy to follow.\n\nOverall, this method appears to be a solid step forward in MCMC for general discrete distributions, well supported with theory and empirical results. \n\n\nWeaknesses:\n\nA key elegance of LB/GWG is the lack of hyper-parameters to tune. MCMC hyper-parameters are notoriously difficult to tune and, at the same time, crucial for good performance. This work adds a path-length proposal distribution as a new hyper-parameter on top of those earlier works. Do the authors have any proposals for how to choose the path-length proposal? I can imagine we'd like to maximize something like the average hop distance. Do the authors have any thoughts on this? If X is set too large for example, then I would assume most proposals would be rejected. There has been much work on this in the context of HMC (to which this method has some spiritual similarities). Even if the authors do not have a concrete approach for choosing this distribution, some discussion of the desiderata of the optimal choice would be appreciated to help guide practitioners in the future. \n\nSection 3.4 should contain a reference to Grathwohl et al. (2021) as this is clearly an extension of ideas presented in that work. As well, theorem 3 of this work is remarkably similar to theorem 1 of that work. This should be noted.\n\nI find the deep EBM (and Ising training) results somewhat misleading. Which PAFS-X do you use? What is the runtime compared to GWG? How many steps are used per iteration in training? If you are using PAFS-X and N steps per iteration then you should compare with GWG using N*E[a(X)] steps since the runtime of PAFS-X is E[a(x)] times the runtime of GWG. How does this change the relative performance of the approaches? When training EBMs one can often get away with simpler MCMC approaches as the model tends to tune itself to the sampler used. This is why SGLD has had much more success in EBMs than more sophisticated samplers such as HMC. If the evaluation was more fair and the performance was still improved I'd be more likely to view this as a straight win over GWG for EBM training. If not, then there would still be some benefit to using GWG since it has no hyper-parameters to tune making it easier for practitioners less familiar with the nuances of tuning MCMC algorithms. \n\nThe authors should add an appendix section on experimental details for the EBM training. How many steps are used per iteration? Which PAFS-X is used? How are the log-likelihood scores evaluated? As it stands the results are not reproducible and I am not convinced of the improvements over GWG presented. Based on the other results in the paper, I believe the results are likely to hold but further evidence should be presented.  \n\nAdd x-axis labels to figures 2 and 4. \n", "summary_of_the_paper": "This paper proposes an MCMC sampling algorithm for discrete distributions. The algorithm can be thought of as an extension of Locally-Balanced proposals from Zanella where the sampler proposes local moves changing one dimension of the input at a time proportional to the likelihood difference of the model in the current state versus each proposed modified state. The key difference here is the authors propose to sample from this proposal distribution multiple times before the accept/reject step. This is somewhat analogous to the difference between Langevin Dynamics and Hamiltonian Monte-Carlo. Ideally, by applying the proposal distribution multiple times before the accept/reject step, the sampler is more likely to escape local modes than the original Locally Balanced algorithm. \n\nThe authors propose a fast variant of this algorithm which is an extension of the recently proposed Gibbs-With-Gradients sampler. This variant replaces the Locally Balanced proposal with the fast approximate version given by Gibbs-With-Gradients. This sampler approximately retains the benefits of the Locally Balanced proposal while being much faster to compute. \n\nThe authors demonstrate the performance of their sampler at sampling from a number of common discrete distributions and at training discrete Energy-Based Models. They find their sampler performs favorably when compared with Locally Balanced proposals, and Gibbs-With-Gradients.  \n\n--------------------- Post rebuttal period ------------------\nI am very glad with the discussion that we've had and I thank the authors for responding to my critiques with new experiments and new results. I think they tell a more complete story of how the method relates to prior work. I believe we have cleared up many issues that other reviewers had with the work. I liked this work initially and I am glad now that other reviewers see it as I did. Having said that, I will keep my recommendation as it was initially. I would be very happy to see this work accepted. \n", "summary_of_the_review": "This work presents a simple extension over prior state-of-the-art discrete MCMC methods which the authors demonstrate leads to improved performance on a wide variety of tasks. The method is simple, easy to implement, and appears correct. The authors provide the necessary theory for understanding the method and provide sufficient evidence for the performance of the method. There are some minor issues with the presentation of this evidence that should be addressed, particularly in the EBM training section. I would support acceptance of this work.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "GZMHOEJoEu9", "reviewer": "Reviewer_Eoj5", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "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": "Locally informed proposals and the use of locally balanced functions were proposed by Zanella (2020). The paper uses these ideas in a new way to efficiently explore large neighborhoods via auxiliary paths.\n\nStrengths:\n1. The method proposed is more efficient (in terms of likelihood evaluations) than naively expanding the neighborhood using a locally-informed proposal.\n2. Using a balanced function the resulting acceptance rate takes a very simple form and can be computed efficiently.\n3. Experiments include many models and several baselines (though I have some comments on the specifics).\n\nWeaknesses:\n1. I find section 3.3 a bit confusing. A few concrete comments: (a) Where does the definition of an \"ideal\" function come from? Theorem 2 states that, for every ideal function f, there is a locally balanced function g asymptotically better than f. Why does this mean asymptotic optimality for locally balanced functions? Am I missing something simple? (b) Is the conditional independence graph (thm 2) defined anywhere? Is it the DAG that defines the factorization of the distribution as in a directed graphical model?\n\n2. The result from theorem 3 (about the use of the Taylor expansion) and its analysis seem problematic. Specifically the part that claims that the use of the approximation does not lead to a large detrimental effect. The theorem lower bounds the efficiency lost by using the approximation. But this lower bound decreases *exponentially* with K and with $U^2$, where U is the length of the path. If U=10 this bound says that using the Taylor approximation leads to a method whose efficiency is at least $T\\times e^{-50}$, where $T$ is the efficiency of the method without the approximation. Claiming that the approximation does not lead to \"losing too much proposal quality\" does not seem justified. It could be the case that this bound is very lose, and thus not informative. Or it could be that there are models for which the inequality becomes an equality, and the use of the Taylor approximation becomes extremely problematic (would it be possible to find a simple model that satisfies this?). Do you have any comments on this? I am not claiming the approximation is necessarily bad, I'm just asking about this specific theoretical result and its relevance.\n\nFinally, I have some comments on the empirical evaluation. First, I'd like to mention that I find the wide variety of models, tasks and baselines very nice. However, I think the section can be improved by showing more organized experiments addressing each component of the proposed method separately. Concretely:\n\nE1. Section 3.2 shows that the method with a path of length L is more efficient than running LB-1 for L steps. It would be nice to have experiments that show the effect of this in practice. This could be done by comparing PAS-L and LB-1 on several models, in terms of likelihood evaluations, not MCMC steps (PAS-L uses L times as many likelihood evaluations per step). Figs 1 and 2 compare these methods in terms of MCMC steps on two simple models. I think the results in terms of likelihood evaluations would be similar. However, these models, while simple and illustrative, may be too tailored to make LB-1 fail. It would be great to see results for the Ising model, RBM, FHMM. \n\nE2. Effect of introducing the Taylor approximation. This could be seen by compare PAS-M vs PAFS-M for a bunch of models, for varying M, in terms of MCMC steps (not likelihood evaluations). While not a fair comparison, this reflects the effect the approximation has, and how much we are losing by using it. (In this setup, we can only expect that using the approximation will hurt performance. This would show how much.) For a practical comparison, results in terms of likelihood evaluations should be shown too (and are currently shown). Also, for the real models, PAS and PAFS are compared only for the Ising model, which the authors state is close to linear (and thus possibly favors the Taylor expansion approximation).\n\nE3. Comparison of the final method with the Taylor approximation against other baselines (GWG-M, LB, Gibbs, etc) in terms of likelihood evaluations.\n\n\nAdditional comments:\n- Should add labels (in x and y axis) in all plots (e.g. x-axis in Fig 4, is it likelihood evaluation? MCMC steps?).\n- Plots in Fig. 3 would be nice to have in same scale for the y-axis.\n- Is there a typo in equation 10? Should it be p instead of n in the summation?", "summary_of_the_paper": "MCMC in discrete spaces using locally informed proposals usually use small neighborhoods (e.g. changing only one component), for efficiency reasons. The paper proposes a method to efficiently explore bigger neighborhoods using locally informed proposals (with a locally balanced function) to sample \"auxiliary paths\". Additionally, following previous work, they propose an efficient version using a linear Taylor expansion, thus saving a significant number of likelihood evaluations.", "summary_of_the_review": "On the positive side, I think the paper provides a nice and natural extension of MCMC methods in discrete spaces based on locally informed proposals. \n\nOn the negative side, I think that some of the theory could be presented more clearly, some of the claims revisited, and the empirical evaluation could be modified a bit to better show the benefits/drawbacks of each of the method's components.\n\n-----\n\nI updated my score after the authors' response.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "AtkYDriAuNu", "reviewer": "Reviewer_y5iJ", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "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": "\n## Weakness\n\n1 - The primary weakness in the paper is that the proof of Theorem 1 is incorrect. Equation (17) can be split into two parts:\n\nA - $\\pi(x)\\sum_{L} \\alpha(L) \\sum_{(\\sigma,L) \\in \\Sigma(\\mathcal{X}, \\mathcal{N}) \\colon \\sigma_0 = x, \\sigma_L = y} \\Big[ \\prod_{l=1}^L Q_0(\\sigma_{l-1},\\sigma_l)\\Big]$\n\nB - $\\min{} \\Big[1, \\frac{\\pi(y)\\prod_{l=1}^L Q_0(\\sigma_l, \\sigma_{l-1})}{\\pi(x)\\prod_{l=1}^L Q_0(\\sigma_{l-1}, \\sigma_{l})}\\Big]$\n\nWhile part A marked above is acurate, part B should technically be  $\\min{} \\Big[      1, \\frac{\\pi(y)\\prod_{l=1}^{L'} Q_0(\\delta_l, \\delta_{l-1})}{\\pi(x)\\prod_{l=1}^{L'} Q_0(\\delta_{l-1}, \\delta_{l})} \\Big]$\n\nwhere $L'$ and $\\delta$ are the length samples and path samples drawn in the algorithm mentioned above Theorem 1.\nAs such Equation 17 doesn't imply Equation 18.\n\nDespite this theoretical inconsistency, it is not surprising that the experimental section shows that the proposed method works better than [GWG]. \nThis is because the method is a Pseudo Marginal MH procedure (where the acceptance probability is a Monte Carlo sample of the true acceptance probability like [Section 4.2, RB]).\nIf so presented, I would be inclined to accept the paper but in the current form, the paper is mathematically incorrect.\n\n2 - The second listed contribution (linearization) is a rehash of [GWG] and is not novel.\n\n3 - There are multiple grammatical errors in the paper examples include:\n  - The first phrase of the abstract should've been `EBMs offer`\n  - Second sentence, second para, page 3 should be `Instead of directly sampling from a large neghborhood path auxiliary sampler sequentially samples new states from a ...`.\n\n4- Minor question: In the thrid paragraph, do you mean a simple random walk? If so what is the graph over which the RW takes place? Did you mean Gibbs sampling?\n\n## Stregths\n\n1 - The experimental section is very strong and the explanation in terms of the toy example is very intuitive.\n\n[GWG] - Will Grathwohl, Kevin Swersky, Milad Hashemi, David Duvenaud, and Chris J Maddison. Oops i took a gradient: Scalable sampling for discrete distributions. arXiv preprint arXiv:2102.04509, 2021.\n[RB] - On Markov chain Monte Carlo methods for tall data\nR Bardenet, A Doucet, C Holmes - The Journal of Machine Learning Research, 2017 https://www.jmlr.org/papers/volume18/15-205/15-205.pdf\n", "summary_of_the_paper": "- Consider the task of MH-MCMC on a EBM defined over discrete state space.\n- [GZ] proposed a Kernel which evaluates the energy function over a Hamming ball surrounding the current state and [GWG] improved the efficiency of the latter by using a Taylor series expansion of the energy function.\n- The authors of the current paper improve the acceptance ratio (Eq 2) and mixing (Thm3) by sampling a multiple transitions (a path) instead of a single neighborhood sample.\n\n[GZ] - Giacomo Zanella. Informed proposals for local mcmc in discrete spaces. Journal of the American Statistical Association, 115(530):852–865, 2020.\n\n[GWG] - Will Grathwohl, Kevin Swersky, Milad Hashemi, David Duvenaud, and Chris J Maddison. Oops i took a gradient: Scalable sampling for discrete distributions. arXiv preprint arXiv:2102.04509, 2021.\n", "summary_of_the_review": "\nThe proof of the main theorem is incorrect and hence I am rejecting the paper.\n\n--------------\nI am updating my scores based on author discussion. Please refer to comments for details.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "2VBplDiZiAF", "reviewer": "Reviewer_EvqS", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "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": "**Contribution and Potential:**\n\nThe work tackles an important open problem for MCMC in discrete spaces, as it enables making large sampling steps while retaining the local benefits of pointwise informed proposals. This has important consequences, for instance when sampling from or training energy-based models. The result in Theorem 2 is to my knowledge a simple but non-trivial extension of the result in (Zanella 2020) and it’s significant in this context, as it provides an initial answer to the question “are balancing functions optimal on large neighborhoods?”. (Indeed, (Zanella 2020) focused on the case of L=1, while the authors consider the extended setting with L>=1.)\n\n**Related Work**\n\nSome recent work must be discussed. For instance, see (Dai 2020), (Jaini 2021) and (Sansone 2021).\n\n**Questions:**\n\n1. Theory. While certainly the result in Theorem 2 holds in the limit of infinite dimension. Is there a limit on the result of Theorem 2, or equivalently an upper bound on the value of L when the dimension of the space is finite? Indeed, as (Zanella 2020) points out, balancing proposals may not be optimal when considering global moves. Similarly, the constant in the Peskun ordering relation (in Theorem 2) should vanish when increasing L. Can the authors elaborate more on that?\n2. Experiments. How is the maximum L determined and how can it be chosen in practice?\n3. Experiments. How is the burn-in time chosen in all experiments?\n4. Related Work. How does the proposed solution relate to existing recent work? \n\n**Minor:**\n\nPage 2 should $\\{0,1\\}^p$ be replaced by $\\{0,1\\}^n$, in order to be consistent with what follows?\n\nPage 7 Figure 10 -> Figure 3\n\nSpacing between text and citations\n\n**References**\n\n(Zanella 2020) Informed Proposals for Local MCMC in Discrete Spaces\n\n(Grathwohl 2021) Oops I Took a Gradient: Scalable Sampling for Discrete Distributions\n\n(Dai 2020) Learning Discrete Energy-based Models via Auxiliary-variable Local Exploration\n\n(Jaini 2021) Sampling in Combinatorial Spaces using SurVAE Flow Augmented MCMC\n\n(Sansone 2021) LSB: Local Self-Balancing MCMC in Discrete Spaces", "summary_of_the_paper": "The paper provides a generalisation of the work of pointwise informed proposals (Zanella 2020), where the acceptance criterion is applied after L sampling steps, instead of each single one.\nThis comes with the advantage that the overall acceptance criterion  is guaranteed to be larger than the corresponding one in (Zanella 2020). The authors also extend the results in (Zanella 2020), by proving under mild conditions that weighting functions based on balancing ones are asymptotically optimal as for sampling efficiency according to Peskun ordering. Lastly, they provide an algorithmic extension reducing the number of target likelihood evaluations by applying the approximation strategy proposed in (Grathwohl 2021). The analysis is corroborated by experiments on several energy based models both for inference and learning tasks, thus showing the superiority of the proposed algorithm and its extension over existing MCMC strategies.\n", "summary_of_the_review": "**Initial recommendation**\n\nThe paper is well written and technically sound. I went through the proofs and confirm their validity.\nThe paper represents an added value to the existing body of knowledge on MCMC, especially in the discrete domain. However, I have some questions for the authors, which might affect the final result. I initially recommend for a weak accept.\n\n**Final recommendation**\n\nThe authors have addressed all concerns during the rebuttal. Therefore, I increase my score.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper provides a novel path auxiliary algorithm for more efficiently exploring discrete state spaces within a Metropolis-Hastings sampler for energy based models. In particular, it essentially replaces the "single site update" by instead proposing an entire path using local information, thus enabling the chain to take larger steps, which can improve acceptance/mixing significantly as they demonstrate. The work is a timely contribution that improves upon exciting recent work. After much discussion among several knowledgeable reviewers and clarifications regarding some details of the main theorem from the authors, there is consensus that the contributions are correct, novel, and likely of impact to the machine learning community. Since the revision period, the empirical evaluations have also been improved and the contributions have methodological novelty as well as promising practical performance.
Decision: Accept (Spotlight)