COMPARISON OF HIT-AND-RUN, SLICE SAMPLER AND RANDOM
WALK METROPOLIS

arXiv:1505.00579v2 [math.PR] 21 Jun 2017

DANIEL RUDOLF AND MARIO ULLRICH
Abstract. Different Markov chains can be used for approximate sampling of a distribution given by an unnormalized density function with respect to the Lebesgue measure.
The hit-and-run, (hybrid) slice sampler and random walk Metropolis algorithm are popular tools to simulate such Markov chains. We develop a general approach to compare the
efficiency of these sampling procedures by the use of a partial ordering of their Markov
operators, the covariance ordering. In particular, we show that the hit-and-run and the
simple slice sampler are more efficient than a hybrid slice sampler based on hit-and-run
which, itself, is more efficient than a (lazy) random walk Metropolis algorithm.

1. Introduction
In many scenarios of Bayesian statistics, statistical physics and other branches of applied
sciences, see [RR02, LV06a, BGJM11], it is of interest to sample on Rd with respect to
a distribution π. In particular, we assume that π is given by an unnormalized density.
More precisely, let K ⊆ Rd be an open, measurable set and ρ : K → (0, ∞) be a positive,
bounded and with respect to the Lebesgue measure an almost everywhere continuous
function. Then, define the probability measure π through ρ by
R
ρ(x) dx
π(A) = R A
ρ(x) dx
K

for all measurable A ⊆ K.
Maybe the most successful approach is the construction of a suitable Markov chain
which can be used to approximate π. The hit-and-run algorithm, the random walk Metropolis, the simple slice sampler and hybrid slice sampler provide such construction
methods. A crucial question is: Which one of these algorithm should be used?
This question is of course related to the speed of convergence of the Markov chain
sampling and any answer depends very much on the imposed assumptions. In general
it is difficult to derive explicit estimates of this speed of convergence. But, it might be
possible to prove that one algorithm is better than another. This motivates the idea of
the comparison of Markov chains.
The first comparison result of Markov chains is due to Peskun [Pes73]. There, a partial
ordering on finite state spaces is invented, where one transition kernel has higher order
than another one if the former dominates the latter off the diagonal. This order was later
Date: June 22, 2017.

2

DANIEL RUDOLF AND MARIO ULLRICH

extended by Tierney [Tie98] to general state spaces. However, for the Markov chains we
have in mind it seems not possible to use this off-diagonal ordering. We consider a partial
ordering on the set of linear operators, see [Kre89, p. 470]. In the context of Markov chains
this ordering is called covariance ordering, see [Mir01, ML09]. Let L2 (π) be the Hilbert
space of functions with finite stationary variance and assume that P1 , P2 : L2 (π) → L2 (π)
are two self-adjoint linear operators. Then we say, P1 ≤ P2 if and only if
hP1 f, f iπ ≥ hP2 f, f iπ ,

f ∈ L2 (π).

Here the inner-product of L2 (π) is given by
Z
hf1 , f2 iπ =
f1 (x)f2 (x) dπ(x),
K

f1 , f2 ∈ L2 (π).

The operator P1 is called positive if P1 ≥ 0, that is, hP1 f, f iπ ≥ 0 for any f ∈ L2 (π).
Reversible transition kernels of given Markov chains induce self-adjoint Markov operators
and we can compare these operators.
Let P1 , P2 be two Markov operators and assume that P1 ≤ P2 . There are a number
of consequences for the corresponding Markov chains: For example, the spectral gap
of P2 is smaller than that of P1 . The spectral gap of a Markov chain is a quantity
which is closely related to the speed of convergence to π, see [Bax05, Rud12]. Another
example is concerning the stationary asymptotic variance of sample averages. For i =
P
(i)
(i)
(i)
1, 2 let Sn (f ) = n1 nk=1 f (Xk ) with Markov chain (Xk )k∈N starting in stationarity
corresponding to Pi and an arbitrary
R function f ∈ L2 (π). These sample averages give
approximations of the mean Ef = K f dπ and, in [Tie98, Mir01] it is observed that, if
P1 ≤ P2 , then
V(f, P2 ) ≤ V(f, P1 ).
(i)

Here V(f, Pi ) := limn→∞ n · E|Sn (f ) − Ef |2 is the stationary asymptotic variance of
(i)
Sn (f ). In other words, the Markov chain of P2 is (asymptotically) more efficient than
that of P1 . For more details to implications of P1 ≤ P2 we refer to Section 2.1.1.
The goal of this article is to compare the hit-and-run algorithm, the (lazy) random
walk Metropolis, the simple slice sampler and a hybrid slice sampler based on hit-and-run
according to this partial ordering. To do so we develop a systematic approach for the
comparison of Markov chains which can be written by a suitable two step procedure.
The intuition behind our approach is the following. It is a simple and well-known observation that, if the self-adjoint and positive Markov operators P1 and P2 satisfy P1 P2 = P2 ,
then P1 ≤ P2 . However, due to the quite restrictive assumptions, this technique can only
be used in very few cases. A useful generalization of this technique, which was invented
in a special case in [Ull14], can be used whenever we have the representation Pi = RPei R∗
for certain operators R and Pei , i = 1, 2. In this case, and if Pe1 Pe2 = Pe2 , one can also
conclude that P1 ≤ P2 , see Lemma 8. Our comparison inequalities therefore follow once
we have established such representations for the Markov chains under consideration.
Let us briefly explain the algorithms. Roughly, a transition of the hit-and-run algorithm
works as follows. Choose randomly a line through the current state and sample according

COMPARISON OF MARKOV CHAINS

3

to π restricted to this line. Thus, instead of sampling with respect to π on K ⊆ Rd hit-andrun only uses sampling with respect to π on 1-dimensional lines through the state space,
which is feasible in many cases 1. In contrast, the simple slice sampler chooses a suitable ddimensional set, a level set of ρ, depending on the current state and then, samples the next
state of the Markov chain uniformly distributed on this level set. Sampling of the uniform
distribution on a d-dimensional set is often not efficiently implementable. This is why this
Markov chain is mostly of theoretical interest. The hybrid slice sampler we are interested
in overcomes this problem by replacing the uniform sampling by one step of a hit-and-run
algorithm on the level set: First, choose a line through the current state uniformly at
random and then generate the next state uniformly distributed on the intersection of that
line with the level set. We call this procedure hybrid slice sampler based on hit-and-run.
Intuitively, it is clear that the simple slice sampler is better than that hybrid one. The
intuition for the comparison of the hit-and-run algorithm and the hybrid slice sampler
based on hit-and-run is not so obvious. Observe that this particular hybrid sampler can
also be interpreted as choosing first a line and, then, performing a simple slice sampling
step according the distribution of π restricted to that line. This observation leads us to
the fact that the hit-and-run algorithm is better. Finally, let us consider the random
walk Metropolis. Assume that we have a proposal density q on Rd and let x ∈ K be the
current state. A transition works as follows: Generate z ∈ Rd according the distribution
determined by q and accept x + z as the next state with probability min{1, ρ(x + z)/ρ(x)}
if x + z ∈ K. Otherwise stay at x. It is well known, see [Hig98, RU13], that the random
walk Metropolis can be interpreted as a certain (hybrid) slice sampling procedure, which
runs a random walk according to q on the level set with uniform limit distribution. We
want to compare the random walk Metropolis and the hybrid slice sampler based on a
hit-and-run. Thus, the question is whether the uniform hit-and-run step is better than
the random walk step on the level set. It turns out that this is indeed the case.
Now let us formulate the main results. To compare the above Markov chains we develop
in Section 3 a general approach which might be of interest on its own. There, in Lemma 8
conditions for two suitably defined Markov operators P1 , P2 are stated which imply that
P1 ≤ P2 . This lemma is the main ingredient for the comparison argument. Its application
leads us to Theorem 9. For the Markov operators M, U, H, S of the (lazy) random walk
Metropolis with rotational invariant proposal q, the hybrid slice sampler based on hitand-run, the hit-and-run algorithm and the simple slice sampler, respectively, we have
M ≤ U ≤ H,
M ≤ U ≤ S.

Thus, the random walk Metropolis is less efficient than the hit-and-run algorithm and
simple slice sampler. The hybrid slice sampler based on hit-and-run we propose lies
concerning efficiency in between.
1For

example, if ρ is a log-concave density, then the restriction of π to a line has also a log-concave
density and one can use different acceptance/rejection methods to sample such a distribution on the line
efficiently. For details we refer to [CR04, Sect. 2.4.2] see also [LV07, Rud13].

4

DANIEL RUDOLF AND MARIO ULLRICH

The paper is organized as follows. In the next section we introduce the notation we use,
comment on the partial ordering, present consequences for the Markov chains and state
the algorithms we study in detail. In Section 3 we invent a new approach how to compare
Markov chains with a specific structure. Section 4 contains the application of the former
developed comparison arguments. Finally, in Section 5, we give some concluding remarks
and discuss open problems.
2. Preliminaries
1/2

Let L2 (π) be the Hilbert space of functions f : K → R with finite norm kf kπ = hf, f iπ ,
where the inner-product of f1 , f2 ∈ L2 (π) is denoted by
Z
hf1 , f2 iπ =
f1 (x)f2 (x) dπ(x).
K

Let P be a transition kernel on K which is reversible with respect to π and let (Xn )n∈N
be a Markov chain with transition kernel P , i.e.
P (x, A) = P(Xk+1 ∈ A | Xk = x)
almost surely for all k ∈ N and A ⊆ K. The corresponding Markov operator, also denoted
by P , is given by
Z
(1)
P f (x) =
f (y)P (x, dy)
K

for f ∈ L2 (π). Obviously, P (x, A) = P 1A (x) for all A ⊆ K, where 1A denotes the
indicator function of A. Note that, by the reversibility, P : L2 (π) → L2 (π) is self-adjoint.
We say that a (Markov) operator P on L2 (π) is positive if hP f, f iπ ≥ 0 for all f ∈ L2 (π).
2.1. On the ordering and consequences. With this notation we define on the set of
Markov operators the following partial ordering. For Markov operators P1 and P2 we
write
P1 ≤ P2
if and only if hP1 f, f iπ ≥ hP2 f, f iπ for all f ∈ L2 (π). In the following let us motivate why
the consideration of this ordering is particularly meaningful for Markov chains.

2.1.1. Consequences for the speed of convergence. There are many ways to measure the
speed of convergence of the distribution of a Markov chain towards its stationary distribution. Probably the most desirable quantity is the total variation distance of νP n and
π, i.e.
kνP n − πktv = sup |νP n (A) − π(A)|,
A⊆K
R
n
n
where νP (A) = K P (x, A) dν(x) is the distribution of the Markov chain with transition
kernel P and initial distribution ν after n steps. However, estimating the total variation
distance is quite delicate and, in practice, it is usually much easier to derive bounds on

COMPARISON OF MARKOV CHAINS

5

it by more analytic quantities, like isoperimetric constants or certain norms of P , see
e.g. [Lov99, LV06b] and [MN07]. Many of these quantities are defined by
cM (P ) = inf h(I − P )f, f iπ
f ∈M

for certain sets of functions M ⊂ L2 (π), where If := f . Hence, the proof of P1 ≤ P2 is
enough to obtain cM (P1 ) ≤ cM (P2 ) for every choice of M. Prominent examples are the
R
• spectral gap: M = {f : kf kπ = 1, K f dπ = 0}
• conductance: M = {f : f = √1A , π(A) ∈ (0, 1/2], A ⊂ K}
π(A)
R
• log-Sobolev constant: M = {f : K f 2 log(f 2 /kf k2π ) dπ = 1}.
There are some more quantities of this form, like the best constant in a Nash inequality,
average and blocking conductance. For details see e.g. [Che05] and [MT06].
In what follows, we will prove P1 ≤ P2 for a couple of Markov operators and, hence,
that P2 is “faster” than P1 in all the above senses.

2.1.2. Consequences for sample averages. The property P1 ≤ P2 has also consequences
for the worst case mean square error and the asymptotic stationary variance of Markov
chain Monte Carlo methods.
(i)
Let (Xk )k∈N be a Markov chain with transition kernel Pi , i = 1, 2, and initial distribution π, and define the Markov chain Monte Carlo method
n

1X
(i)
f (Xk ),
i = 1, 2,
n k=1
R
which gives an approximation to Ef := K f dπ for f ∈ L2 (π). By virtue of spectral
theoretic arguments2 one can show that
Sn(i) (f ) =

sup E Sn(2) (f ) − Ef

kf kπ ≤1

2

≤ sup E Sn(1) (f ) − Ef

2

.

kf kπ ≤1

Another interesting consequence is stated in [ML09, Theorem 6]. There it is proven
that P1 ≤ P2 if and only if
V(f, P2 ) ≤ V(f, P1 ),
(i)

for each individual f ∈ L2 (π), where V(f, Pi ) := limn→∞ n · E|Sn (f ) − Ef |2 is the
asymptotic stationary variance (which can also be considered as asymptotic mean square
error).
2.2. The algorithms. We present the different algorithms we consider in detail and
provide relevant literature.
2One

argues with Corollary 3.27 and Lemma 2.12 from [Rud12], as well as the fact that the spectral
gap of P1 is smaller than that of P2 .

6

DANIEL RUDOLF AND MARIO ULLRICH

2.2.1. Hit-and-run algorithm. The hit-and-run algorithm proposed by Smith [Smi84] is
well studied in different settings, see [BRS93, DF97, KS98, Lov99, LV06b, KSZ11, Rud13].
Informally it samples at each step on a randomly chosen 1-dimensional line with respect
to the corresponding conditional distribution. Let Sd−1 be the Euclidean unit sphere in
Rd . For x ∈ K and θ ∈ Sd−1 we define
L(x, θ) = {x + sθ ∈ K | s ∈ R}

as the chord through x in direction θ. A transition from x to y of the hit-and-run algorithm
works as follows: Generate a set L(x, θ) by choosing θ with the uniform distribution on
the sphere and, then, choose y ∈ L(x, θ) according to the distribution determined by ρ
conditioned on the chord L(x, θ). A transition of hit-and-run is given by Algorithm 1.
Algorithm 1. (Hit-and-run)
Transition from current state x to next state y:
(1) Sample θ ∼ Uniform(Sd−1 ).
(2) Sample y ∼ Hθ (x, ·), where
R
1 (z)ρ(z) dz
L(x,θ) A
R
Hθ (x, A) =
.
ρ(z) dz
L(x,θ)

Note that the integral in Hθ is over a 1-dimensional subset of Rd and the integration
is with respect to the 1-dimensional Lebesgue measure. For x ∈ K and A ⊆ K the
transition kernel, say H, of the hit-and-run algorithm is determined by
Z
H(x, A) =
Hθ (x, A) dσ(θ),
Sd−1

where σ = Uniform(Sd−1 ) denotes the uniform distribution on the sphere. It is well known
that this transition kernel is reversible with respect to π, see for example [BRS93].
An important special case of the hit-and-run algorithm above is given if the density
is an indicator function, for example ρ = 1K . In this case, the hit-and-run algorithm is
reversible with respect to the uniform distribution on K. Thus, under weak regularity
assumptions, the uniform distribution is the (unique) stationary distribution, see [BRS93].
We call this special case uniform hit-and-run. Let us mention that we use the uniform
hit-and-run in the next section.
2.2.2. Simple and hybrid slice sampling. Slice sampling belongs to the class of auxiliary
variable algorithms that are defined by a Markov chain on an extended state space, see
[Hig98, RR99, MMR01, MT02, RR02, Nea03] and the references therein.
We consider the simple slice sampler and a hybrid slice sampler based on hit-and-run.
A single transition of the simple slice sampler is presented in Algorithm 2.
Algorithm 2. (Simple slice sampler)
Transition from current state x to next state y.
(1) Sample t ∼ Uniform(0, ρ(x)).

COMPARISON OF MARKOV CHAINS

7

(2) Sample y ∼ Uniform(K(t)), where

K(t) = {x ∈ K | ρ(x) > t}

is the level set of ρ determined by t.

The transition kernel, say S, corresponding to Algorithm 2 is
Z ρ(x)
vold (A ∩ K(t))
1
dt,
S(x, A) =
ρ(x) 0
vold (K(t))
where vold denotes the d-dimensional Lebesgue measure. The simple slice sampler exhibits
quite robust convergence properties, see [RR99, MT02]. However, a crucial drawback is
that the second step is difficult to implement. Because of this we consider the following
hybrid slice sampler based on hit-and-run. The idea is to replace the second step of the
simple slice sampler by a Markov chain transition according to the uniform hit-and-run
algorithm in K(t), see Algorithm 3.
Algorithm 3. (Hybrid slice sampler based on hit-and-run)
Transition from current state x to next state y.
(1) Sample t ∼ Uniform(0, ρ(x)) and θ ∼ Uniform(Sd−1 ) independently.
(2) Sample y ∼ Uniform (Lt (x, θ)), where
Lt (x, θ) = {x + rθ ∈ K(t) | r ∈ R}

is the chord through x in direction θ restricted to K(t).
The transition kernel, say U, of this hybrid slice sampler is given by
Z ρ(x) Z
1
|Lt (x, θ) ∩ A|
U(x, A) =
dσ(θ) dt
ρ(x) 0
|Lt (x, θ)|
Sd−1

where |·| here denotes the 1-dimensional Lebesgue measure. This modification is tempting
since the uniform hit-and-run algorithm is easily implementable, at least when the level
sets are convex, i.e. ρ is quasi-concave. Further, the simple slice sampler and the hybrid
slice sampler are reversible with respect to π, see for example in [LR14].
2.2.3. Random walk Metropolis. The random walk Metropolis in Rd provides an easy to
implement method for Markov chain sampling. Further, it is well studied and different
convergence results are known, see for example [RT96, MT96, JH00].
To guarantee that a certain operator is positive we consider a lazy version of a random
walk Metropolis. Further, we assume in the following that q : Rd → [0, ∞) is a rotational
invariant probability density on Rd , i.e. q(rθ1 ) = q(rθ2 ) for θ1 , θ2 ∈ Sd−1 . The rotational
invariance guarantees that q is symmetric. Let us provide two examples which satisfy the
rotational invariance.
Example 4. Let x ∈ Rd , δ > 0 and Bd be the Euclidean unit ball with κd := vold (Bd ).
Then we can set q(x) := 1Bd (δ −1 x)/(δ d κd ) and the corresponding random walk Metropolis
is known as δ-ball walk.

8

DANIEL RUDOLF AND MARIO ULLRICH

Example 5. Again, let x ∈ Rd and set q(x) := exp(− |x|2 /2)/(2π)d/2 . The corresponding
random walk Metropolis is known as Gaussian random walk.
In Example 4 the proposal q depends on a parameter δ. In connection to this we want to
mention that in recent years optimal scaling results concerning a parameter of the proposal
of the random walk Metropolis attracted a lot of attention, see [RGG97, RR01, NRY12].
Now a single transition of the (lazy) random walk Metropolis with proposal q is described in Algorithm 6.
Algorithm 6. (Random walk Metropolis)
Transition from current state x to next state y with proposal q.
(1) Sample z ∼ q and u1 , u2 ∼ Uniform (0, 1) independently.
(2) If x + z ∈ K, u1 ≤ 1/2 and u2 < min{1, ρ(x + z)/ρ(x)} then accept the proposal,
and set y := x + z, else reject the proposal and set y := x.
To simplify the notation we define the acceptance probability of a proposed state y =
x + z as
n ρ(y) o
α(x, y) = min 1,
ρ(x)
for x, y ∈ K and α(x, y) = 0 otherwise. Then, the transition kernel is given by
Z
1
M(x, A) =
1A (y) α(x, y) q(y − x) dy,
2 K
for A ⊆ K with x ∈
/ A and M(x, {x}) = 1 − M(x, K \ {x}).
3. Auxiliary variable Markov chains
We develop a systematic approach how to compare Markov chains which can be described by a suitable two-step procedure. The idea of a comparison of this type was
developed in [Ull12, Ull14] in a specific setting.
For this let A be an arbitrary (index) set equipped with a measure λ. For example
A = Sd−1 and λ can be the spherical measure or A = [0, ∞) and λ is just the Lebesgue
measure. Further, we have a function s : K × A → [0, ∞) which satisfies:
• for all x ∈ K we have that s(x, ·) is a probability density function according to λ;
• for all a ∈ A we have that s(·, a) is integrable according to π.
For every a ∈ A we assume that we have an equivalence relation ∼a and a transition
kernel Pa on (K, B(K)), such that Pa (x, A) = 0 for each x ∈ K and A ⊆ K \ [x]a , where
[x]a = {y ∈ K : x ∼a y} is the equivalence class of K by ∼a to which x belongs.
We consider Markov chains in K for which a single transition, starting from x ∈ K,
can be written as the following procedure:
1) Sample a ∈ A according to the distribution with density s(x, ·).
2) Generate the next state with respect to Pa (x, ·).

COMPARISON OF MARKOV CHAINS

9

That is, the transition kernel P is given by
Z

(2)
P (x, A) =
Pa x, A s(x, a) dλ(a)
A

for A ⊆ K.

Remark 7. Since the support of the measure Pa (x, ·) is contained in [x]a and [x]a = [y]a
if y ∈ [x]a , we can interpret the Markov chain with transition kernel Pa on K and initial
state x as a Markov chain on [x]a .
Clearly, for every Markov chain a transition can be written in this form. (Simply, set
Pa (x, ·) := P (x, ·) and take [x]a ≡ K for arbitrary A and λ.) However, in some cases
it is possible to represent the transitions of two different Markov chains in the form (2)
with the same equivalence relations ∼a and s(x, a), so that the corresponding “inner”
transition kernels are much easier to analyze. We show that a suitable relation of these
“inner” kernels is enough to compare the original Markov chains.
To state the result it is convenient to define for almost all a ∈ A (concerning λ) a
probability measure πa on K induced by s(x, a), i.e.
R
s(x, a) dπ(x)
πa (A) = R A
,
A ⊆ K.
s(x, a) dπ(x)
K
(1)

(2)

Lemma 8. Assume that for all a ∈ A there are transition kernels Pa , Pa
that
(1)
(2)
(i) Pa and Pa are reversible with respect to πa ;
(1)
(ii) Pa is positive on L2 (πa );
(1) (2)
(2)
(iii) Pa Pa = Pa .
Then, for the operators P1 , P2 : L2 (π) → L2 (π) defined by
Z
Pi f (x) =
Pa(i) f (x) s(x, a) dλ(a), i ∈ {1, 2}

on K such

A

holds P1 ≤ P2 .
Proof. First, we show that Pi can be written as a product of suitable operators. For this
let µ be a probability measure on K × A that is given by
Z Z
(3)
µ(B) :=
1B (x, a) s(x, a) dλ(a) dπ(x)
K

A

for B ⊆ K × A, and let L2 (µ) be the Hilbert space of functions g : K × A → R with finite
1/2
norm kgkµ = hg, giµ , where the inner-product of g, h ∈ L2 (µ) is defined by
Z
hg, hiµ =
g(x, a) h(x, a) dµ(x, a)
K×A
Z Z
=
g(x, a) h(x, a) s(x, a) dλ(a) dπ(x).
K

A

10

DANIEL RUDOLF AND MARIO ULLRICH

Further, let R : L2 (µ) → L2 (π) be the operator that is given by
Z
(4)
Rg(x) =
g(x, a) s(x, a) dλ(a),
g ∈ L2 (µ).
A

Since
hf, Rgiπ =

Z Z
K

f (x) g(x, a) s(x, a) dλ(a) dπ(x),
A

we obtain that the adjoint operator R∗ : L2 (π) → L2 (µ) satisfies
R∗ f (x, a) = f (x),

(5)

f ∈ L2 (π),

for a ∈ A and x ∈ K. Finally, with the operator Pei : L2 (µ) → L2 (µ) that is given by
Z
(6)
Pei g(x, a) =
g(y, a) Pa(i)(x, dy),
g ∈ L2 (µ),
K

we obtain

Pi = RPei R∗ ,

(7)

i ∈ {1, 2}.

We want to prove different properties of Pe1 and Pe2 . For this it is useful to write
ga (x) = g(x, a) for g ∈ L2 (µ) and fixed a ∈ A. Clearly, ga ∈ L2 (πa ) for almost every
(i)
a ∈ A (w.r.t. λ), such that Pa ga is well-defined. For i ∈ {1, 2} we show that Pei is
(i)
self-adjoint on L2 (µ). Note that by the assumptions we know that Pa is self-adjoint on
L2 (πa ). With
Z
Ca =
s(x, a) dπ(x)
K

we have for g, h ∈ L2 (µ) that
Z Z Z
e
hPi g, hiµ =
g(y, a) Pa(i)(x, dy) h(x, a) s(x, a) dλ(a) dπ(x)
K A K
Z
Z Z
=
Ca
g(y, a) Pa(i)(x, dy) h(x, a) dπa (x) dλ(a)
K K
Z
ZA
(i)
Ca hga , Pa(i) ha (·, a)iπa dλ(a)
=
Ca hPa ga , ha iπa dλ(a) =
A

A

= hg, Peihiµ .

(1)

By the same line of arguments one can show that Pe1 preserves the positivity of Pa , i.e.
Pe1 is positive on L2 (µ). Further, for g ∈ L2 (µ) we obtain
Pe1 Pe2 g(x, a) = Pa(1) Pa(2) ga (x) = Pa(2) ga (x) = Pe2 g(x, a).

By the self-adjointness of Pe1 and Pe2 we also have Pe2 Pe1 = Pe2 .
Now we gathered all tools together to prove the assertion. By the positivity, we know
that Pe1 has a unique positive square root N, i.e. N 2 = Pe1 . It is well known that N

COMPARISON OF MARKOV CHAINS

11

commutes with every operator that commutes with Pe1 , see e.g. [Kre89, Theorem 9.4-2],
in particular N Pe2 = Pe2 N. We obtain
hP2 f, f iπ = hRPe2 R∗ f, f iπ = hPe2 R∗ f, R∗ f iµ
= hPe1 Pe2 R∗ f, R∗ f iµ = hPe2 NR∗ f, NR∗ f iµ

≤ hNR∗ f, NR∗ f iµ = hPe1 R∗ f, R∗ f iµ = hP1 f, f iπ ,

where the inequality comes from kPe2 k ≤ 1, which is true since Pe2 is a Markov operator.

Thus, for each of the intended comparisons, say between P1 and P2 , we have to
• find representations
Z
Pi f (x) =
Pa(i) f (x)s(x, a) dλ(a), i ∈ {1, 2};
A

(i)

• check the assumptions (i)-(iii) of Lemma 8, i.e. reversibility of Pa with respect to
(1)
(1) (2)
(2)
πa , positivity of Pa and Pa Pa = Pa .
4. Main result
In this section we apply Lemma 8 to prove the following theorem.
Theorem 9. Let M, U, H, S be the Markov operators of the (lazy) random walk Metropolis
with rotational invariant proposal q, the hybrid slice sampler based on hit-and-run, the hitand-run algorithm and the simple slice sampler. Then
M ≤ U ≤ H,
M ≤ U ≤ S.
Before proving the different inequalities in the statement we start with recalling some
notion and state a useful lemma. Recall that Sd−1 denotes the Euclidean unit sphere in
Rd and σ is the uniform distribution on Sd−1 . For θ ∈ Sd−1 , x ∈ Rd and t ∈ [0, ∞) let
K(t) = {x ∈ K | ρ(x) > t}

be the level set of ρ determined by t and

Lt (x, θ) = {x + θr ∈ K(t) | r ∈ R}

be the chord in K(t) through x in direction θ. Moreover, define L(x, θ) := L0 (x, θ) and,
for a set A ⊆ K = K(0), let
Πθ (A) = {y ∈ Rd | y⊥θ, L(y, θ) ∩ A 6= ∅}

be the orthogonal projection of A to the hyperplane that is orthogonal to θ. We obtain
the following useful results by an application of the Theorem of Fubini and the integral
transformation to polar coordinates.

12

DANIEL RUDOLF AND MARIO ULLRICH

Lemma 10. Let t ≥ 0 and θ ∈ Sd−1 . For Lebesgue integrable f : K(t) → R we have
Z
Z
Z
(8)
f (x) dx =
f (y) dy dx,
K(t)

Πθ (K(t))

Lt (x,θ)

and for fixed x ∈ Rd holds
Z
Z
Z
dκd
f (y) |x − y|d−1 dy dσ(θ)
(9)
f (y) dy =
2 Sd−1 Lt (x,θ)
K(t)
where κd = vold (Bd ) denotes the volume of the d-dimensional Euclidean unit ball.
Note that in both identities the inner integral on the ride-hand-side is over a 1-dimensional
subset of Rd and the integration is with respect to the 1-dimensional Lebesgue measure.
The different inequalities in Theorem 9 will be proven in the following sections. There
we use the notation of Section 3.
4.1. Hit-and-run vs. hybrid slice sampler: Let A := Sd−1 and λ := σ. Further, for
(x, a) ∈ K × A set s(x, a) := 1 such as [x]a := L0 (x, a). This implies πa = π for all a ∈ A.
By
R
ρ(y)1A (y) dy
[x]a
R
Ha (x, A) =
,
A ⊆ K,
ρ(y)
dy
[x]a
we can write the Markov operator H : L2 (π) → L2 (π) of the hit-and-run algorithm as
Z Z
Hf (x) =
f (y) Ha(x, dy) dλ(a).
A

Let
1
Ua (x, A) =
ρ(x)

K

Z

ρ(x)

0

Z

Lt (x,a)

1A (y)
dy dt
|Lt (x, a)|

and observe that Lt (x, a) = [x]a ∩ K(t). With this we can write the Markov operator
U : L2 (π) → L2 (π) of the hybrid slice sampler based on hit-and-run as
Z Z
Uf (x) =
f (y) Ua(x, dy) dλ(a).
A

K

Thus we have a common representation.
Lemma 8:

Now we check the assumptions (i)-(iii) of

(i) To the reversibility of Ha and Ua with respect to πa :
We only show reversibility of Ua , since reversibility of Ha follows by the same line of
arguments. For A, B ⊆ K it is enough to prove
Z
Z
(10)
Ua (x, B)ρ(x) dx =
Ua (x, A)ρ(x) dx.
A

B

COMPARISON OF MARKOV CHAINS

13

By the application of (8) with s = 0, the equality 1[0,ρ(y)) (t) = 1K(t) (y) and the fact
that y ∈ [x]a implies Lt (y, a) = Lt (x, a) we obtain
Z
Z
Z
1A (y)Ua (y, B)ρ(y) dy dx
Ua (x, B)ρ(x) dx =
A

=

Z

Πa (K)

=

Z

Πa (K)

Πa (K)

Z

1A (y)

[x]a
∞

Z

0

Z

0

[x]a

ρ(y)

Z

Lt (y,a)

1B (z)
dz dt dy dx
|Lt (y, a)|

1
|Lt (x, a) ∩ A| |Lt (x, a) ∩ B| dt dx.
|Lt (x, a)|

Observe that in the right-hand-side A and B are interchangeable which proves (10).
(ii) To the positivity of Ua :
R
By similar arguments as in the proof of (10) we obtain with c = K ρ(x) dx that
2
Z
Z
Z ∞
1
1
f (y) dy
dt dx ≥ 0.
hUa f, f iπ =
c Πa (K) 0 |Lt (x, a)|
Lt (x,a)
(iii) To Ua Ha = Ha :
From the definition of the Markov kernels Ha it is obvious that y ∈ [x]a implies
Ha (y, A) = Ha (x, A) for all A ⊆ K. This implies
Z ρ(x) Z
1
dy
Ua Ha (x, A) =
Ha (y, A)
dt = Ha (x, A)
ρ(x) 0
|Lt (x, a)|
Lt (x,a)

and hence Ua Ha = Ha .
Thus, as a direct consequence of Lemma 8 we obtain U ≤ H.

4.2. Simple vs. hybrid slice sampler: Here, we derive another representation of the
hybrid slice sampler adapted to simple slice sampling.
For this let A := [0, ∞) and λ be the 1-dimensional Lebesgue measure. Further, for
(x, a) ∈ K × A set s(x, a) := 1[0,ρ(x)) (a)/ρ(x) such as [x]a := K(a). Observe that
πa (A) =
For any x ∈ K let

vold (A ∩ K(a))
.
vold (K(a))

Sa (x, A) = πa (A),

A ⊆ K.

Clearly, the Markov operator S : L2 (π) → L2 (π) of the simple slice sampler can be written
as
Z Z
Sf (x) =

f (y) Sa(x, dy) s(x, a) dλ(a).

A

K

Note that, with these notations, we have La (x, θ) = [x]a ∩ L(x, θ) and let
Z
Z
1A (y)
Ua (x, A) =
dy dσ(θ).
Sd−1 La (x,θ) |La (x, θ)|

14

DANIEL RUDOLF AND MARIO ULLRICH

With this we can write the Markov operator U : L2 (π) → L2 (π) of the hybrid slice sampler
based on hit-and-run as
Z Z
Uf (x) =
f (y) Ua(x, dy)s(x, a) dλ(a).
A

K

Again, we have a common representation and it remains to check the assumptions (i)-(iii)
of Lemma 8:
(i) To the reversibility of Sa and Ua with respect to πa :
Since Sa (x, A) = πa (A) reversibility of Sa is obvious. Observe that πa is the uniform
distribution on [x]a = K(a) and Ua (x, ·) performs a uniform hit-and-run step on
K(a) which is known to be reversible, see for example [Rud12, Lemma 4.10].
(ii) To the positivity of Ua :
With c = vold (K(a)), by the application of (8) and by the fact that for y ∈ La (x, θ)
follows La (x, θ) = La (y, θ) holds
Z
Z
Z
1
f (z)f (x)
hUa f , f iπa =
dz dσ(θ) dx
c K(a) Sd−1 La (x,θ) |La (x, θ)|
Z
Z
Z
Z
f (y)f (z)
1
dz dy dx dσ(θ)
=
c Sd−1 Πθ (K(a)) La (x,θ) La (y,θ) |La (y, θ)|
Z
2
Z
Z
1
1
=
f (y) dy
dx dσ(θ) ≥ 0.
c Sd−1 Πθ (K(a)) |La (x, θ)|
La (x,θ)
(iii) To Ua Sa = Sa :
This follows immediately from Sa = πa and the reversibility of Ua with respect to
πa .
Thus, as a direct consequence of Lemma 8 we obtain U ≤ S.
4.3. Hybrid slice sampler vs. Metropolis: Again we need a suitable representation
for the hybrid slice sampler based on hit-and-run.
Here, let A := Sd−1 × [0, ∞) and λ be the product measure of σ and the 1-dimensional
Lebesgue measure. Further, for x ∈ K and (a1 , a2 ) ∈ A, set s(x, a1 , a2 ) = 1[0,ρ(x)] (a2 )/ρ(x)
such as [x](a1 ,a2 ) := La2 (x, a1 ). Observe that
π(a1 ,a2 ) (A) =

vold (A ∩ K(a2 ))
,
vold (K(a2 ))

is the uniform distribution in K(a2 ). By
Z
U(a1 ,a2 ) (x, A) =

[x](a1 ,a2 )

1A (y)

A ⊆ K,
dy
[x](a1 ,a2 )

we have a representation of the Markov operator U : L2 (π) → L2 (π) of the hybrid slice
sampler by
Z Z
Uf (x) =
f (y)Ua (x, dy)s(x, a) dλ(a)
A

K

COMPARISON OF MARKOV CHAINS

15

with a ∈ A. Now we have to represent the random walk Metropolis in the same fashion.
For this let us define
dκd
η(x, y) =
q(y − x) |y − x|d−1 , x, y ∈ Rd
2
and for A ⊆ K let


Z
Z
1
1
1A (y) η(x, y) dy + 1A (x) 1 −
η(x, y) dy .
Ma (x, A) =
2 [x]a
2 [x]a
Here
it is essential that q is rotational invariant, namely this property assures that
R
η(x, y) dy ≤ 1. Note that Ma is again a lazy transition kernel, since Ma (x, {x}) ≥ 1/2.
[x]a
For x 6∈ A by (9) we have
Z
1
M(x, A) =
1A (y) α(x, y) q(y − x) dy
2 K
Z ρ(x) Z
1
1[0,ρ(y)] (t)1A (y) q(y − x) dy dt
=
2ρ(x) 0
K
Z ρ(x) Z
1
=
2
1A (y) q(y − x) dy dt
2ρ(x) 0
K(t)
Z ρ(x) Z
Z
dκd
=
1A (y)q(x − y) |x − y|d−1 dy dt dσ(θ)
4ρ(x) Sd−1 0
Lt (x,θ)
Z
=
Ma (x, A)s(x, a) dλ(a).
A

Thus, the transition kernel of the random walk Metropolis has the desired representation.
It remains to check the conditions (i)-(iii) of Lemma 8:
(i) To the reversibility of Ma and Ua with respect to πa :
This follows again by a suitable application of (8). Let a = (a1 , a2 ). Due to its simple
form, reversibility of Ua is obvious. It is enough to show for disjoint A, B ⊆ K that
Z
Z
Ma (x, B) dx =
Ma (x, A) dx.
A∩K(a2 )

B∩K(a2 )

Since La2 (x, a1 ) = La2 (y, a1 ) if y ∈ La2 (x, a1 ) we have
Z
Z
1
1B (y)η(x, y) dy dx
2 A∩K(a2 ) [x]a
Z
Z
1
=
1A (x)1B (y)η(x, y) dy dx
2 K(a2 ) La2 (x,a1 )
Z
Z
Z
1
=
1A (y)1B (z)η(y, z) dz dy dx
2 Πa1 (K(a2 )) La2 (x,a1 ) La2 (y,a1 )
Z
Z
Z
1
=
1A (y)1B (z)η(y, z) dz dy dx.
2 Πa1 (K(a2 )) La2 (x,a1 ) La2 (x,a1 )

16

DANIEL RUDOLF AND MARIO ULLRICH

By the symmetry of q, i.e. q(y − z) = q(z − y) follows η(y, z) = η(z, y). This leads
to the reversibility.
(ii) To the positivity of Ma :
By the fact that Ma is a lazy transition kernel we have positivity.
(iii) To Ma Ua = Ua :
By the fact that y ∈ [x]a implies [x]a = [y]a and hence Ua (x, A) = Ua (y, A) for all
A ⊆ K we obtain
Ma Ua (x, A)


Z
Z
1
1
=
Ua (y, A) η(x, y) dy + Ua (x, A) 1 −
η(x, y) dy
2 [x]a
2 [x]a
= Ua (x, A).
Thus, as a direct consequence of Lemma 8 we obtain M ≤ U.
5. Concluding remarks
In this article we have presented a technique to compare the efficiency of Markov chains
of a specific type. Using this technique we provide two comparison hierarchies according
to a partial ordering of Markov operators of four prominent Markov chains for sampling
general distributions in Rd . The comparison with respect to the partial ordering leads
to comparison results according to different criterions, for example the spectral gap, the
conductance or the log-Sobolev constant, cf. Section 2.1.1.
Let us mention here that the computational cost for the simulation of each individual
Markov chain is not taken into account. There seems to be a trade-off between efficiency
and computional cost which should be further investigated. We leave this open for future
work.
Finally, there are three open problems related to the considered Markov chains.
First, what is about the relation of the hit-and-run and the simple slice sampler? It is
easy to see that there cannot be a general result as in the other cases. For this consider
two examples:
1. If π is the uniform distribution on K, then one step of the simple slice sampler is
enough to sample π, while for the hit-and-run algorithm it is not (as long as d > 1).
Hence, in this situation, hit-and-run is worse.
2. Let d = 1. Then, the hit-and-run algorithm samples π in one step, regardless of π.
Hence, in this situation hit-and-run is better.
It seems to be interesting to find cases where hit-and-run is better. This is because, we
guess that in general the computational cost of the simple slice sampler is (if it is at all
implementable) much higher.
The second open problem concerns reverse inequalities. That is, if a Markov chain is
better than another, how much better is it? This seems to be a delicate question and
we can answer it only for toy examples. The techniques that were used in [Ull14] in a
discrete setting do not seem to work here.

COMPARISON OF MARKOV CHAINS

17

The last problem we want to mention: Is there a similar hierarchy for the mixing
time? That is, the number of steps that are needed to make the total variation distance
“small”, cf. Section 2.1.1. Certainly, the answer to this question additionally depends on
the choice of the initial distribution. But the (quite analytical) techniques of this paper,
see also [Pes73, LR14, Ull14, AV14], do not seem to be suitable for this purpose. One
interesting approach in this direction for Markov chains on discrete state spaces is given
in [FK13].
References
[AV14]

Ch. Andrieu and M. Vihola, Establishing some order amongst exact approximations of
MCMCs, Preprint, Available at http://arxiv.org/abs/1404.6909 (2014).
[Bax05]
P. Baxendale, Renewal theory and computable convergence rates for geometrically ergodic
Markov chains, Ann. Appl. Probab. 15 (2005), no. 1B, 700–738.
[BGJM11] S. Brooks, A. Gelman, G. Jones, and X. Meng (eds.), Handbook of markov chain monte carlo,
Chapman and Hall, 2011.
[BRS93] C. Bélisle, H. Romeijn, and R. Smith, Hit-and-run algorithms for generating multivariate
distributions, Math. Oper. Res. 18 (1993), no. 2, 255–266.
[Che05]
M. Chen, Eigenvalues, inequalities, and ergodic theory, Probability and its Applications (New
York), Springer-Verlag London Ltd., London, 2005.
[CR04]
G. Casella and C. Robert, Monte carlo statistical methods, second edn., Springer Texts in
Statistics (New York), Springer-Verlag, 2004.
[DF97]
P. Diaconis and D. Freedman, On the hit and run process, Tech. rep. 497, Department of
Statistics, University of California, Berkeley (1997), 1–16.
[FK13]
J. Fill and J. Kahn, Comparison inequalities and fastest mixing markov chains, Ann. App.
Probab. 23 (2013), 1778–1816.
[Hig98]
D. Higdon, Auxiliary variable methods for markov chain monte carlo with applications, J.
Amer. Statist. Assoc. 93 (1998), no. 442, 585–595.
[JH00]
S. Jarner and E. Hansen, Geometric ergodicity of Metropolis algorithms, Stochastic Process.
Appl. 85 (2000), no. 2, 341–361.
[Kre89]
E. Kreyszig, Introductory functional analysis with applications, John Wiley & Sons, Inc., 1989.
[KS98]
D. Kaufman and R. Smith, Direction choice for accelerated convergence in hit-and-run sampling, Oper. Res. 46 (1998), no. 1, 84–95.
[KSZ11] S. Kiatsupaibul, R. Smith, and Z. Zabinsky, An analysis of a variation of hit-and-run for
uniform sampling from general regions, ACM Trans. Model. Comput. Simul. 21 (2011), no. 3,
16:1–16:11.
[Lov99]
L. Lovász, Hit-and-run mixes fast, Math. Program. 86 (1999), no. 3, Ser. A, 443–461.
[LR14]
K. Latuszyński and D. Rudolf, Convergence of hybrid slice sampling via spectral gap, Arxiv
(2014), submitted.
[LV06a]
L. Lovász and S. Vempala, Fast algorithms for logconcave functions: sampling, rounding, integration and optimization, Proceedings of the 47th Annual IEEE Symposium on Foundations
of Computer Science (Washington, DC, USA), FOCS ’06, IEEE Computer Society, 2006,
pp. 57–68.
, Hit-and-run from a corner, SIAM J. Comput. 35 (2006), no. 4, 985–1005.
[LV06b]
[LV07]
, The geometry of logconcave functions and sampling algorithms, Random Structures
Algorithms 30 (2007), no. 3, 307–358.
[Mir01]
A. Mira, Ordering and improving the performance of Monte Carlo Markov chains, Statist.
Sci. 16 (2001), no. 4, 340–350.

18

DANIEL RUDOLF AND MARIO ULLRICH

[ML09]

A. Mira and F. Leisen, Covariance ordering for discrete and continuous time Markov chains,
Statistica Sinica 19 (2009), 651–666.
[MMR01] A. Mira, J. Møller, and G. Roberts, Perfect slice samplers, J. Roy. Statist. Soc. Ser. B 63
(2001), no. 3, 593–606.
[MN07]
P. Mathé and E. Novak, Simple Monte Carlo and the Metropolis algorithm, J. Complexity 23
(2007), no. 4-6, 673–696.
[MT96]
K. Mengersen and R. Tweedie, Rates of convergence of the Hastings and Metropolis algorithms,
Ann. Statist. 24 (1996), no. 1, 101–121.
[MT02]
A. Mira and L. Tierney, Efficiency and convergence properties of slice samplers, Scand. J.
Statist. 29 (2002), no. 1, 1–12.
[MT06]
R. Montenegro and P. Tetali, Mathematical aspects of mixing times in Markov chains, Found.
Trends Theor. Comput. Sci. 1 (2006), no. 3, x+121. MR 2341319 (2009g:68100)
[Nea03]
R. Neal, Slice sampling, Ann. Statist. 31 (2003), no. 3, 705–767, With discussions and a
rejoinder by the author.
[NRY12] P. Neal, G. Roberts, and W. Yuen, Optimal scaling of random walk Metropolis algorithms
with discontinuous target densities, Ann. Appl. Probab. 22 (2012), no. 5, 1880–1927.
[Pes73]
P. Peskun, Optimum Monte-Carlo sampling using Markov chains, Biometrika 60 (1973), 607–
612.
[RGG97] G. Roberts, A. Gelman, and W. Gilks, Weak convergence and optimal scaling of random walk
Metropolis algorithms, Ann. Appl. Probab. 7 (1997), no. 1, 110–120.
[RR99]
G. Roberts and J. Rosenthal, Convergence of slice sampler Markov chains, J. R. Stat. Soc.
Ser. B Stat. Methodol. 61 (1999), no. 3, 643–660.
, Optimal scaling for various Metropolis-Hastings algorithms, Statist. Sci. 16 (2001),
[RR01]
no. 4, 351–367.
, The polar slice sampler, Stoch. Models 18 (2002), no. 2, 257–280.
[RR02]
[RT96]
G. Roberts and R. Tweedie, Geometric convergence and central limit theorems for multidimensional Hastings and Metropolis algorithms, Biometrika 83 (1996), no. 1, 95–110.
[RU13]
D. Rudolf and M. Ullrich, Positivity of hit-and-run and related algorithms, Electron. Commun.
Probab. 18 (2013), 1–8.
[Rud12]
D. Rudolf, Explicit error bounds for Markov chain Monte Carlo, Dissertationes Math. 485
(2012), 93 pp.
, Hit-and-run for numerical integration, J. Dick, F. Y. Kuo, G. Peters, I. H. Sloan
[Rud13]
(Eds.), Monte Carlo and Quasi-Monte Carlo Methods 2012, Springer Proceedings in Mathematics & Statistics 65 (2013), 597–612.
[Smi84]
R. Smith, Efficient Monte Carlo procedures for generating points uniformly distributed over
bounded regions, Oper. Res. 32 (1984), no. 6, 1296–1308.
[Tie98]
L. Tierney, A note on Metropolis-Hastings kernels for general state spaces, Ann. Appl. Probab.
8 (1998), no. 1, 1–9.
[Ull12]
M. Ullrich, Rapid mixing of Swendsen-Wang and single-bond dynamics in two dimensions,
ArXiv e-prints (2012).
, Rapid mixing of Swendsen-Wang dynamics in two dimensions, Dissertationes Math.
[Ull14]
502 (2014), 64 pp.
Institute for Mathematical Stochastics, Universität Göttingen
E-mail address: daniel.rudolf@uni-goettingen.de
URL: http://stochastik.math.uni-goettingen.de/~drudolf/
Institut für Analysis, Johannes Kepler Universität Linz
E-mail address: mario.ullrich@jku.at

