# Fundamentals of Partial Rejection Sampling

Mark Jerrum  School of Mathematical Sciences, Queen Mary, University of
London, Mile End Road, London E1 4NS

###### Abstract.

Partial Rejection Sampling is an algorithmic approach to obtaining a perfect
sample from a specified distribution. The objects to be sampled are assumed to
be represented by a number of random variables. In contrast to classical
rejection sampling, in which all variables are resampled until a feasible
solution is found, partial rejection sampling aims at greater efficiency by
resampling only a subset of variables that ‘go wrong’. Partial rejection
sampling is closely related to Moser and Tardos’ algorithmic version of the
Lovász Local Lemma, but with the additional requirement that a specified
output distribution should be met. This article provides a largely self-
contained account of the basic form of the algorithm and its analysis.

This work was supported by EPSRC grant EP/S016694/1, “Sampling in hereditary
classes”.

##  1\. The setting

The aim of this expository article is to provide a uniform treatment of a
particular approach to sampling combinatorial structures. The method is a
development of classical rejection sampling. Suppose $\Phi(\mathbf{X})$ is a
predicate (Boolean function) depending on random variables
$\mathbf{X}=(X_{1},\ldots,X_{n})$ coming from a product distribution. We would
like to obtain a sample from the conditional distribution of $\mathbf{X}$
given that $\Phi(\mathbf{X})$ holds. Classical rejection sampling repeatedly
generates realisations of $\mathbf{X}$ from the product distribution until one
that satisfies $\Phi$ is found, and then outputs that. In many situations this
approach is very inefficient, as satisfying assignments to $\Phi$ may occur
with exponentially small probability (in $n$). The idea in Partial Rejection
Sampling (PRS) is to identify small subsets of the variable that violate
$\Phi$ (in some sense), and resample just those variables. Clearly, the choice
of which variables to resample has to be done with great care if the output
distribution is to remain the correct one. PRS was identified as a general
algorithmic technique by Guo, Jerrum and Liu [14].

We focus on examples of PRS inspired by algorithmic proofs of the Lovász Local
Lemma (LLL). There is a substantial literature on this topic to assist us, but
it is concerned only with problem of constructing _some_ satisfying assignment
to $\Phi$. Partial rejection sampling adds the novel requirement that the
output should be uniform (or, more generally, from the desired distribution)
on satisfying assignments. This additional requirement adds a new challenge.

In the study of the LLL, the class of ‘extremal’ instances receives particular
attention. The extremal instances are particularly suited to PRS and we treat
them first. It is unlikely that anything in this section of the article is
conceptually new. However, the treatment of PRS in the extremal regime
involves some particularly beautiful combinatorial ideas, and it is seems a
good time to bring this material together in one place, with a consistent
approach and notation.

After that, we investigate to what extent the conditions defining extremal
instances can be relaxed. The viewpoint taken in this section is novel to a
certain extent. Care has been taken to set out the conditions under which PRS
continues to function in the non-extremal setting, in the hope that it will
help in discovering new applications. More than usual attention is paid to the
flexibility in the order in which variables can be resampled.

The scope of this article is limited to versions of partial rejection sampling
that stay close to the spirit of the algorithmic LLL pioneered by Moser and
Tardos [25]. In particular, the number of random variables is finite and all
constraints on them are ‘hard’. Informally, we restrict attention to a
‘combinatorial’ setting, which excludes important applications to spin systems
in statistical physics. We finish with a few pointers to work that goes beyond
the framework presented here.

It should be noted that PRS is not the only approach to perfect sampling.
Perhaps the best known and most extensively studied is ‘Coupling From The
Past’ (CFTP), which was pioneered by Propp and Wilson [27]. Other approaches
include Fill and Huber’s ‘randomness recycler’ [8].

##  2\. Partial rejection sampling in the extremal setting

Suppose $\Phi=\varphi_{1}\wedge\varphi_{2}\wedge\cdots\wedge\varphi_{m}$ is
formula on variables $\mathbf{X}=(X_{1},X_{2},\ldots,X_{n})$. Each clause
$\varphi_{k}$ depends on a certain tuple $(X_{i_{1}},\ldots,X_{i_{a_{k}}})$ of
variables, where $a_{k}$ is the _arity_ of $\varphi_{k}$. We refer to the
tuple ${\sf Scp}(\varphi_{k})=(i_{1},\ldots,i_{a_{k}})$ of indices as the
_scope_ of the clause $\varphi_{k}$. We imagine that the arities $\\{a_{k}\\}$
are small in relation to $n$, and possibly uniformly bounded. By assuming
$i_{1}<i_{2}<\cdots<i_{a_{k}}$ we can consider the scope ${\sf
Scp}(\varphi_{k})$ to be a $a_{k}$-tuple or a set of cardinality $a_{k}$,
according to context. For a set $S\subseteq\\{1,\ldots,n\\}$ of indices we
write $X_{S}=\\{X_{i}:i\in S\\}$. Then to emphasise the dependence on the
variables we can write

$\Phi(\mathbf{X})=\varphi_{1}(X_{{\sf
Scp}(\varphi_{1})})\wedge\cdots\wedge\varphi_{m}(X_{{\sf Scp}(\varphi_{m})}).$

The formula $\Phi$ is said to be an _extremal instance_ if, for all $1\leq
k<\ell\leq m$ satisfying ${\sf Scp}(\varphi_{k})\cap{\sf
Scp}(\varphi_{\ell})\not=\emptyset$ it is the case that $\varphi_{k}(X_{{\sf
Scp}(\varphi_{k})})\vee\varphi_{\ell}(X_{{\sf Scp}(\varphi_{\ell})})$ is a
tautology. In other words, any two clauses that are probabilistically
dependent cannot both be false. In this section we consider only extremal
instances, as these can be dealt with using the basic form of partial
rejection sampling.

Now suppose that variable $X_{i}$ takes values from a set $D_{i}$. Equip
$D_{i}$ with a probability distribution and call the probability space
$\mathcal{D}_{i}$. We are interested in sampling a realisation of the random
variable $\mathbf{X}$ from the product distribution
$\mathcal{D}=\mathcal{D}_{1}\times\mathcal{D}_{2}\times\cdots\times\mathcal{D}_{n}$
conditioned on $\Phi(\mathbf{X})$ holding. Denote this desired distribution by
$\mathcal{D}_{\Phi}$. Partial Rejection Sampling (PRS) is a simple algorithm
for accomplishing this task in the context of extremal instances. It runs as
follows.

Algorithm 1 Partial Rejection Sampling

$\mathrm{PRS}(\Phi,\mathcal{D})$ // $\Phi$ is a formula on variable set
$\mathbf{X}$

Sample $\mathbf{X}$ from the product distribution
$\mathcal{D}=\mathcal{D}_{1}\times\cdots\times\mathcal{D}_{n}$

while $\neg\Phi(\mathbf{X})$ do

Choose any clause $\varphi_{k}$ with $\neg\varphi_{k}(X_{{\sf
Scp}(\varphi_{k})})$

Resample all variables in ${\sf Scp}(\varphi_{k})$

end while

In the resampling step, we are naturally using the product distribution
$\prod_{i\in{\sf Scp}(\varphi_{k})}\mathcal{D}_{i}$.

The algorithm PRS was first introduced by Moser and Tardos [25] in the context
of an algorithmic proof of the Lovász Local Lemma. Its application to sampling
from naturally specified distributions was proposed by Guo, Jerrum and Liu
[14], who analysed its correctness and efficiency. Although their
investigation seems to be the first attempt to treat PRS as a general
technique, several specific examples had previously appeared in the
literature, for example Propp and Wilson’s ‘cycle popping’ algorithm for
spanning trees [27] and Cohn, Pemantle and Propp’s ‘sink-popping’ algorithm
for sink-free orientations [4].

######  Remark.

In the combinatorial community, the Moser and Tardos algorithm would be viewed
as gradually eliminating the set of ‘bad events’ until none are left. In the
area of constraint satisfaction, the goal is to simultaneously satisfy a
collection of constraints. It is important to keep in mind that, of these
diametrically opposing conventions, we use the latter here.

In classical rejection sampling we would resample the whole of $\mathbf{X}$ on
each iteration. In contrast, PRS resamples only a subset of offending
variables. We cannot expect the correctness of the algorithm to survive such
extreme corner cutting. Indeed, for general formulas $\Phi$, the call
$\mathrm{PRS}(\Phi,\mathcal{D})$ does not produce a sample from the
distribution $\mathcal{D}_{\Phi}$. Surprisingly, PRS _does_ achieve the
desired distribution for extremal instances.

######  Theorem 1.

Suppose $\Phi$ is a satisfiable extremal instance. Then
$\mathrm{PRS}(\Phi,\mathcal{D})$ terminates with probability 1. On
termination, $\mathbf{X}$ is a realisation of a random variable from the
distribution $\mathcal{D}_{\Phi}$.

To analyse the algorithm, we need to introduce time explicitly. A _resampling
table_ is a semi-infinite matrix $(\omega_{i,j}:1\leq i\leq n\text{ and
}j\in\mathbb{N})$. Each entry $\omega_{i,j}$ in the table is an independent
sample from the distribution $\mathcal{D}_{i}$. Fixing $i$, the sequence
$\omega_{i,0},\omega_{i,1},\omega_{i,2},\ldots$ specifies the sequence of
values taken by the random variable $X_{i}$ during the execution of the
algorithm. Introducing a superscript to indicate the time $t$ (each iteration
of the loop takes one time unit), we write $X_{i}^{t}=\omega_{i,j(i,t)}$. If
$X_{i}$ is resampled during iteration $t$ then $j(i,t)=j(i,t-1)+1$, otherwise
$j(i,t)=j(i,t-1)$. Initially, $j(i,0)=0$ for all $1\leq i\leq n$. (By
convention, iteration $t$ occupies the interval between time $t-1$ and time
$t$.) At any time $t$, the _frontier_ of the resampling table is
$F(t)=(j(1,t),j(2,t),\ldots,j(n,t))$. See Figure 1.

$\vdots$ | $\vdots$ | $\vdots$ | $\vdots$ | $\vdots$
---|---|---|---|---
$\omega_{1,4}$ | $\omega_{2,4}$ | $\omega_{3,4}$ | $\omega_{4,4}$ | $\omega_{5,4}$
$\omega_{1,3}$ | $\omega_{2,3}$ | $\omega_{3,3}$ | $\omega_{4,3}$ | $\omega_{5,3}$
$\omega_{1,2}$ | $\omega_{2,2}$ | $\omega_{3,2}$ | $\omega_{4,2}$ | $\omega_{5,2}$
$\omega_{1,1}$ | $\omega_{2,1}$ | $\omega_{3,1}$ | $\omega_{4,1}$ | $\omega_{5,1}$
$\omega_{1,0}$ | $\omega_{2,0}$ | $\omega_{3,0}$ | $\omega_{4,0}$ | $\omega_{5,0}$
$X_{1}$ | $X_{2}$ | $X_{3}$ | $X_{4}$ | $X_{5}$
Figure 1. A resampling table

As time progresses, we record the actions of the algorithm in the form of a
partition of the portion of the resampling table that lies behind the
frontier, namely $(\omega_{i,j}:1\leq i\leq n\text{ and }0\leq j<j(i,t))$. The
partition builds as the frontier advances. In iteration $t$, the variables in
the scope of some clause $\varphi_{k}$ are resampled. The locations
$\\{(i,j(i,t-1)):i\in{\sf Scp}(\varphi_{k}))\\}$ that were on the frontier now
lie behind it; this set of locations now forms a new block of the partition.
We call these blocks the _resampling blocks_. At time $t$, the frontier
together with the partition into resampling blocks forms a _transcript_ of the
run of the algorithm up to time $t$.

By way of example, consider the formula

$\Phi(\mathbf{X})=(X_{1}\vee X_{2})\wedge(\neg X_{1}\vee X_{3}\vee\neg
X_{4})\wedge(\neg X_{2}\vee\neg X_{3}\vee X_{5})\wedge(X_{4}\vee\neg X_{5})$

on variables $\mathbf{X}=(X_{1},X_{2},X_{3},X_{4},X_{5})$. (The formula $\Phi$
encodes sink-free orientations of a certain 4-vertex graph, a point we shall
return to later.) Thus $\varphi_{2}=\neg X_{1}\vee X_{3}\vee\neg X_{4}$ and
${\sf Scp}(\varphi_{2})=\\{1,3,4\\}$, and similarly for the other clauses. A
particular realisation of the resampling table that leads to termination of
the algorithm PRS, together with its associated transcript, are depicted in
Figure 2. In the pictorial representation of the transcript, the values in the
resampling table are spread out along the columns so that each resampling
block of the transcript occupies a single row. The rectangle at the top
denotes the _final frontier_ , i.e., the frontier at termination. Initially,
$X_{4}=0$ and $X_{5}=1$, which violates clause $\varphi_{5}$. Accordingly,
variables $X_{4}$ and $X_{5}$ are resampled, and $\\{(4,0),(5,0)\\}$ becomes
the first resampling block of the transcript. The value of $X_{4}$ switches
from 0 to 1, and this causes clause $\varphi_{2}$ to be violated, since now
$X_{1}=X_{4}=1$ and $X_{3}=0$. So $X_{1}$, $X_{3}$ and $X_{4}$ are resampled
and $\\{(1,0),(3,0),(4,1)\\}$ becomes the next resampling block of the
transcript. Eventually, $\mathbf{X}=(0,0,1,1,1)$ , which satisfies $\Phi$, and
the algorithm halts.

Suppose we run the algorithm PRS twice, using different non-deterministic
choices (of which clauses to resample), until termination. A priori, it might
be imagined that the two runs would in general have different transcripts, but
this is not the case, as we shall see in Lemma 2. Some intuition can be gained
from Figure 2. At time 4, $\mathbf{X}=(0,0,0,0,1)$, and hence clauses
$\varphi_{1}$ and $\varphi_{4}$ are both violated. We can resample either
$\\{X_{1},X_{2}\\}$ first or $\\{X_{4},X_{5}\\}$, but either way we end up
with the _same_ transcript. In this context, it is crucial that ${\sf
Scp}(\varphi_{1})\cap{\sf Scp}(\varphi_{4})=\emptyset$, but in an extremal
instance, this condition is guaranteed.

$\vdots$ |  |  | $\vdots$ |
---|---|---|---|---
0 | $\vdots$ | $\vdots$ | 1 | $\vdots$
1 | 0 | 1 | 1 | 1
0 | 0 | 0 | 0 | 1
0 | 1 | 1 | 1 | 0
1 | 0 | 0 | 0 | 1
$X_{1}$ | $X_{2}$ | $X_{3}$ | $X_{4}$ | $X_{5}$

01101001100001101 $t=7$ $t=6$ $t=4,5$ $t=3$ $t=2$ $t=1$ $t=0$
$X_{1}$$X_{2}$$X_{3}$$X_{4}$$X_{5}$00111

Figure 2. A realisation of a resampling table, and the corresponding
transcript

######  Lemma 2.

Fix a resampling table. Suppose that for some sequence of non-deterministic
choices, $\mathrm{PRS}(\Phi,\mathcal{D})$ terminates with a certain
transcript. Then for any other sequence of choices, the algorithm will
terminate with the same transcript.

To prove this lemma, we use a version of Newman’s Lemma that is particularly
convenient in this application. An (abstract) rewriting system is simply a set
$\mathcal{T}$ of ‘positions’ together with a binary ‘rewriting’ relation $\to$
on $T$. For positions $t,s\in\mathcal{T}$, the relation $t\to s$ indicates
that it is possible to go from $t$ to $s$ in one move. A position $t$ from
which no valid move $t\to s$ is possible is said to be _terminal_. A sequence
of moves ending at a terminal state is said to be _terminating_. Following
Eriksson [5], we say that the rewriting system $(\mathcal{T},\to)$ has the
_polygon property_ if, given any position $t\in\mathcal{T}$ and two moves
$t\to s$ and $t\to s^{\prime}$, either (a) there are two sequences $s=s_{0}\to
s_{1}\to\cdots\to s_{\ell}=t^{*}$ and $s^{\prime}=s^{\prime}_{0}\to
s^{\prime}_{1}\to\cdots\to s^{\prime}_{\ell}=t^{*}$ of the same length $\ell$
that end at the same position $t^{*}$, or (b) there are two infinite sequences
of moves starting from $s$ and $s^{\prime}$. A rewriting system is said to
have the _strong convergence property_ if, for any starting position $t$ from
which there exists a sequence of moves terminating at some position $t^{*}$,
it is the case that every sequence of moves starting from $t$ will lead to
$t^{*}$, and in the same number of moves. Eriksson [5, Thm 2.1] showed the
following.

######  Lemma 3.

A rewriting system has the strong convergence property iff it has the polygon
property.

###### Proof of Lemma 2.

Fix a resampling table. View the collection of all possible transcripts as an
abstract rewriting system by introducing a binary relation $\to$ on
transcripts. The meaning of $t\to s$ is that $s$ can follow $t$ in one
iteration of $\mathrm{PRS}$. This rewriting system has the diamond property,
namely if $t\to s$ and $t\to s^{\prime}$ then there exists $t^{*}$ such that
$s\to t^{*}$ and $s^{\prime}\to t^{*}$. (This diamond property is clearly
stronger than the polygon property.) For suppose $t\to s$ is a result of
resampling the variables in scope ${\sf Scp}(\varphi_{k})$, and $t\to
s^{\prime}$ the result of resampling ${\sf Scp}(\varphi_{\ell})$. Since the
instance $\Phi$ is extremal we know that ${\sf Scp}(\varphi_{k})\cap{\sf
Scp}(\varphi_{\ell})=\emptyset$. Thus, we can resample whichever scope was not
resampled in the first step, to get to a common transcript $t^{*}$. The result
now follows from Lemma 3. ∎

###### Proof of Theorem 1.

Fix a particular satisfying assignment
$\mathbf{X}=\mathbf{b}=(b_{1},\ldots,b_{n})$ to $\Phi$. At any point in the
execution of the algorithm, the following fortuitous sequence of events may
occur over the next $n$ iterations: each time a variable $X_{i}$ is resampled,
it is assigned the value $b_{i}$. On each iteration, $\mathbf{X}$ approaches
closer to $\mathbf{b}$ in Hamming distance. Thus, the algorithm will terminate
in the next $n$ iterations. Since this fortuitous sequence of events occurs
with probability bounded away from 0, the running time of the algorithm PRS is
stochastically dominated by an exponential random variable with finite mean.
So the algorithm terminates with probability 1.

Fix a resampling table $T$, and run $\mathrm{PRS}$ on $T$ to obtain a
transcript. Since the algorithm has terminated, we know that the frontier
contains a satisfying assignment. Create a new resampling table $T^{\prime}$
by replacing the values in the frontier by some other satisfying assignment.
Now run the algorithm again on $T^{\prime}$ with the same nondeterministic
choices of scopes to resample. Note that this is always possible: whenever the
algorithm running on $T$ resamples ${\sf Scp}(\varphi_{k})$ at time $t$ it is
because $\varphi_{k}(\omega_{i,j(i,t)}:i\in{\sf Scp}(\varphi_{k}))$ is false.
None of the resampled variables are in the final frontier, since no variables
beyond the final frontier are ever inspected. So the clause $\varphi_{k}$ is
also false when the algorithm is run on table $T^{\prime}$, and it is valid
step to resample ${\sf Scp}(\varphi_{k})$. Finally, on the same iteration that
the algorithm terminates when run on table $T$, it will also terminate on
$T^{\prime}$. The same transcript (i.e., frontier $F(t)$ together with the
partition of the table behind the frontier) arises from running the algorithm
on $T^{\prime}$ as the one that arose from the run on $T$.

By Lemma 2 any sequence of non-deterministic choices made by the algorithm on
table $T^{\prime}$ leads to the same transcript. Summarising, the final
transcript does not depend on the nondeterministic choices made by the
algorithm, and is also unchanged if one satisfying assignment is substituted
for another in the final frontier. Thus, conditioned on the transcript, each
satisfying assignment $\mathbf{X}=(b_{1},\ldots,b_{n})$ of $\Phi$ occurs with
probability proportional to
$\mathcal{D}_{1}(b_{1})\mathcal{D}_{2}(b_{2})\cdots\mathcal{D}_{n}(b_{n})$.
So, at termination, $\mathbf{X}$ is distributed as $\mathcal{D}_{\Phi}$. ∎

There is a remarkably simple (though not simple to derive) formula for the
expected runtime of algorithm PRS, which we now present. Given $\Phi$, define
$\Gamma=\Gamma(\Phi)$ to be the _dependency graph_ with vertex set $[m]$
(where vertex $k$ corresponds to clause $\varphi_{k}$) and edge relation
$\sim$ defined by $k\sim\ell$ iff ${\sf Scp}(\varphi_{k})\cap{\sf
Scp}(\varphi_{\ell})\not=\emptyset$. Let
$\Sigma=\\{s_{1},s_{2},\ldots,s_{m}\\}$ be an alphabet of $m$ symbols. If
$k\sim\ell$ then symbols $s_{k}$ and $s_{\ell}$ do not commute; otherwise,
$s_{k}$ and $s_{\ell}$ do commute, i.e., $s_{k}s_{\ell}=s_{\ell}s_{k}$. Denote
by $\mathcal{R}_{\Gamma}$ the set of commutation relations:

$\mathcal{R}_{\Gamma}=\\{s_{k}s_{\ell}=s_{\ell}s_{k}:k\not\sim\ell\\}.$

The set $\Sigma^{*}/\mathcal{R}_{\Gamma}$ of _traces_ over $\Sigma$ is the set
of all words over the alphabet $\Sigma$ quotiented by the commutation
relations $\mathcal{R}_{\Gamma}$. So a trace can be thought of as word over
$\Sigma$ where we regard two words as indistinguishable if one can be obtained
from the other by transposing adjacent commuting symbols.

Cartier and Foata [3] derived the generating function for traces. Introduce
indeterminates $z_{1},\ldots,z_{m}$ corresponding to the $m$ clauses in
$\Phi$, and define

$P_{\Gamma}(z_{1},\ldots,z_{m})=\sum_{I\in\mathcal{I}(\Gamma)}(-1)^{|I|}z_{I},$

where $z_{I}=\prod_{k\in I}z_{k}$ and $\mathcal{I}(\Gamma)$ is the set of all
independent sets in $\Gamma$. Note that the polynomial $P_{\Gamma}$ is the
generating function of independent sets in $\Gamma$, with terms signed
according to parity. Cartier and Foata showed that
$P_{\Gamma}(z_{1},\ldots,z_{m})^{-1}$ is the generating function for traces
$\Sigma^{*}/\mathcal{R}$. (The derivation can also be found, e.g., in Knuth
[20, Thm F] and Viennot [30, Prop. 5.1].)

######  Lemma 4.

The generating function $T_{\Gamma}$ for traces $\Sigma^{*}/\mathcal{R}$ is
given by $T_{\Gamma}(z_{1},\ldots,z_{m})=P_{\Gamma}(z_{1},\ldots,z_{m})^{-1}$.
To be precise, $P_{\Gamma}T_{\Gamma}=1$ in the algebra of traces.

In particular, ignoring the non-commutativity of variables for a moment, the
coefficient of $z_{1}^{e_{1}}z_{2}^{e_{2}}\cdots z_{m}^{e_{m}}$ in
$T_{\Gamma}(z_{1},\ldots,z_{m})$ is the number of traces in which symbol
$s_{1}$ occurs $e_{1}$ times, $s_{2}$ occurs $e_{2}$ times, etc.

The motivation for introducing traces is that they are in perfect
correspondence with transcripts, where the $m$ symbols correspond to the $m$
possible kinds of resampling blocks; specifically, symbol $s_{k}$ corresponds
to a block arising from resampling ${\sf Scp}(\varphi_{k})$. Let
$w=s_{i_{1}}s_{i_{2}}\ldots s_{i_{t}}$ be any word in $\Sigma^{*}$. Consider
the transcript that results if the algorithm PRS performs block resamplings in
the order ${\sf Scp}(\varphi_{i_{1}}),{\sf Scp}(\varphi_{i_{2}}),\ldots,{\sf
Scp}(\varphi_{i_{t}})$. Now let
$w^{\prime}=s_{i_{1}^{\prime}}s_{i_{2}^{\prime}}\ldots s_{i_{t}^{\prime}}$ be
any word in $\Sigma^{*}$ that is equivalent to $w$ under the commutation
relations $\mathcal{R}$. It is not difficult to see that the same transcript
results from the sequence of block resamplings ${\sf
Scp}(\varphi_{i_{1}^{\prime}}),{\sf Scp}(\varphi_{i_{2}^{\prime}}),\ldots,{\sf
Scp}(\varphi_{i_{1}^{\prime}})$. (Transposing the order of two adjacent
commuting symbols transposes the order in which two blocks are resampled;
however, those blocks have no variables in common, so there is no change in
the transcript.) Conversely, if words $w$ and $w^{\prime}$ lead to the same
transcript then they must be equivalent under commutativity. (Suppose
$s_{i_{1}}\not=s_{i_{1}^{\prime}}$. Let $s_{i_{h}^{\prime}}$ be the first
occurrence of the symbol $s_{i_{1}}$ in $w^{\prime}$. The first $h-1$
resamplings prompted by $w^{\prime}$ did not disturb the variables in ${\sf
Scp}(\Phi_{i_{1}})$. Therefore, $s_{i_{h}^{\prime}}$ commutes with all earlier
symbols in $w^{\prime}$ and can be ‘bubbled’ into first place. The remaining
symbols can be brought into alignment inductively.) Transcripts are exactly
the _empilements [des pièces]_ or ‘heaps of pieces’ of Viennot [30], who gives
a beautiful pictorial explanation of the correspondence between _empilements_
(and hence transcripts) and traces. See also Knuth [20, §7.2.2.2].

The correspondence between traces and transcripts can be appreciated
pictorially in Figure 2. The depiction of the transcript is based on Viennot’s
_empilements_. Knuth invites us to think of each symbol as a piece in Tetris
that appears from above and descends until further progress is obstructed. The
word $w=s_{4}s_{2}s_{1}s_{3}s_{1}s_{4}s_{2}$ specifies an order for the
arriving pieces that leads to the transcript on the right of the figure. The
word $w^{\prime}=s_{4}s_{2}s_{1}s_{3}s_{4}s_{1}s_{2}$ leads to the same
transcript, since $s_{1}$ and $s_{4}$ commute. In contrast, the word
$s_{4}s_{2}s_{3}s_{1}s_{1}s_{4}s_{2}$ results in a different transcript, as
$s_{1}$ and $s_{3}$ do not commute: the pieces corresponding to symbols
$s_{1}$ and $s_{3}$ cannot pass each other. The equivalence class
$\\{w,w^{\prime}\\}$ is a trace, since the only adjacent commuting pair of
symbols is $s_{1}s_{4}$. Traces, _empilements_ and transcripts are different
views of the same concept.

Before analysing the runtime of algorithm PRS, let us observe that it is
remarkably easy to compute the probability of observing a particular
transcript such as the one in Figure 2. Recall that each value in the
resampling table is the result of an independent toss of a fair coin. At time
$t=0$, we have that $\omega_{4,0}=0$ and $\omega_{5,0}=1$, an event that
occurs with probability $\frac{1}{4}$. (The only way for
$\varphi_{4}(X_{4},X_{5})$ to be false is for $X_{4}$ to be 0, and $X_{5}$ to
be 1.) At time $t=1$ we know that $\omega_{1,0}=1$, $\omega_{3,0}=0$ and
$\omega_{4,1}=1$, an event with probability $\frac{1}{8}$, and so on for times
$t=2,3,4,5,6$. All these events are independent, and the probability that they
all occur is $2^{-17}$. Finally, the frontier must contain a satisfying
assignment; there are 10 satisfying assignments out of a total of 32, so the
probability of observing the transcript depicted is $10\times 2^{-22}$.

For $k\in[m]$, let
$p_{k}=\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\neg\varphi_{k})$ denote
the probability that $\varphi$ is false in the product distribution, and
extend this notation to a set of clauses $S\subseteq[m]$ by letting
$p_{S}=\prod_{k\in S}p_{k}$. Then define

$q_{S}=\sum_{I\in\mathcal{I}(\Gamma):I\supseteq S}(-1)^{|I\setminus S|}p_{I}.$

Note that $q_{S}=0$ if $S\notin\mathcal{I}(\Gamma)$. Note also that

(1) $\displaystyle P_{\Gamma}(p_{1},\ldots,p_{m})$
$\displaystyle=\sum_{I\in\mathcal{I}(\Gamma)}(-1)^{|I|}p_{I}=q_{\emptyset}$
$\displaystyle p_{k}P_{\Gamma-N[k]}(p_{1},\ldots,p_{m})$
$\displaystyle=p_{k}\sum_{I\in\mathcal{I}(\Gamma-N[k])}(-1)^{|I|}p_{I}$
$\displaystyle=\sum_{I\in\mathcal{I}:I\ni k}(-1)^{|I|-1}p_{I}$ (2)
$\displaystyle=q_{\\{k\\}},$

where $\Gamma-N[k]$ denotes the graph obtained from the dependency graph
$\Gamma$ by removing the closed neighbourhood of $k$ and incident edges. (The
_closed neighbourhood_ $N[k]$ of $k$ is the set containing vertex $k$ and all
its neighbours.)

In the case of extremal instances, the quantity $q_{\emptyset}$ has a simple
probabilistic interptetation. By the principle of inclusion-exclusion,

$\displaystyle\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\Phi)$
$\displaystyle=\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}\bigg{(}\bigwedge_{k\in[m]}\varphi_{k}\bigg{)}$
$\displaystyle=\sum_{S\subseteq[m]}(-1)^{|S|}\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}\bigg{(}\bigwedge_{k\in
S}\neg\varphi_{k}\bigg{)}$ (3)
$\displaystyle=\sum_{I\in\mathcal{I}(\Gamma)}(-1)^{|I|}\prod_{k\in
I}\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\neg\varphi_{k})$ (see below)
$\displaystyle=\sum_{I\in\mathcal{I}(\Gamma)}(-1)^{|I|}p_{I}$ (4)
$\displaystyle=q_{\emptyset}.$

Equality (3) uses two facts: (a) when $S$ is not an independent set the
corresponding term is zero, and (b) for any independent set $I$, the events
$\\{\neg\varphi_{k}:k\in I\\}$ are probabilistically independent. Note, in
particular, that $q_{\emptyset}>0$ when $\Phi$ is satisfiable.

In the next theorem, the first sampling of the variables in some scope is
regarded as a _re_ sampling, even though there was no previous one.

######  Theorem 5.

Suppose $\Phi$ is a satisfiable extremal instance. Then the expected number of
resamplings of the scope of $\varphi_{k}$ during a run of
$\mathrm{PRS}(\Phi,\mathcal{D})$ is $q_{\\{k\\}}/q_{\emptyset}$.

###### Proof.

As noted earlier, the generating function for transcripts is
$P_{\Gamma}(z_{1},\ldots,z_{m})^{-1}$. We claim that the generating function
for transcripts weighted according to probability of occurrence is
$q_{\emptyset}P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})^{-1}$. In other words,
the probability of observing a transcript with $e_{k}$ resamplings of scope
${\sf Scp}(\varphi_{k})$, for $1\leq k\leq m$, is the coefficient of
$z_{1}^{e_{1}}z_{2}^{e_{2}}\cdots z_{m}^{e_{m}}$ in
$q_{\emptyset}P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})^{-1}$. To see this, fix
a transcript with $e_{k}$ resamplings of scope ${\sf Scp}(\varphi_{k})$, for
$1\leq k\leq m$, and consider the probability that a random resampling table
will generate that transcript. The frontier must contain a satisfying
assignment, which happens with probability $q_{\emptyset}$, by (4). Each block
corresponding to a clause $\varphi_{k}$ must contain an assignment making
$\varphi_{k}$ false, which happens with probability $p_{k}$. All these
probabilities are independent, so the overall probability of observing the
transcript is $q_{\emptyset}\prod_{1\leq k\leq m}p_{k}^{e_{k}}$. The claim
follows. Note that we have used that the fact that if it is _possible_ for a
certain transcript to arise from a given resampling table it _will_ do so.

Note that $1-P_{\Gamma}(p_{1},\ldots,p_{m})=1-q_{\emptyset}\in[0,1)$, and so
the power series expansion

$P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})^{-1}=\frac{1}{1-(1-P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m}))}=\sum_{i=0}^{\infty}\big{(}1-P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})\big{)}^{i}$

converges in an open neighbourhood of the point $z_{1}=\cdots=z_{m}=1$. The
expected number of resamplings of the scope of $\varphi_{k}$ is given by

$\displaystyle q_{\emptyset}\frac{\partial}{\partial
z_{k}}P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})^{-1}\bigg{|}_{z_{1}=\cdots=z_{m}=1}$
$\displaystyle\qquad\hbox{}=-q_{\emptyset}P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})^{-2}\frac{\partial}{\partial
z_{k}}P_{\Gamma}(p_{1}z_{1},\ldots,p_{m}z_{m})\bigg{|}_{z_{1}=\cdots=z_{m}=1}$
$\displaystyle\qquad\hbox{}=q_{\emptyset}P_{\Gamma}(p_{1},\ldots,p_{m})^{-2}p_{k}P_{\Gamma-N[k]}(p_{1},\ldots,p_{m}).$

We use here the fact that $P_{\Gamma}$ is multilinear, so differentiating with
respect to $x_{k}$ eliminates terms corresponding to independent sets that do
not include $k$. Using identities (1) and (2), we see that the expected number
of times ${\sf Scp}(\varphi_{k})$ is resampled is $q_{\\{k\\}}/q_{\emptyset}$.
∎

We can recast the above theorem in a simple, easy to use form.

######  Corollary 6.

The expected number of iterations of Algorithm PRS on input
$(\Phi,\mathcal{D})$ is

$\mathop{\mathbb{{}E}}\nolimits(\textup{\\#iterations})=\frac{\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\textup{Exactly
one clause in $\Phi$ is
true})}{\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\textup{$\Phi$ is
true})}.$

###### Proof.

Generalising the inclusion-exclusion argument used earlier, and assuming
$S\in\mathcal{I}(\Gamma)$, we have

$\displaystyle\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}\bigg{(}\bigwedge_{k\in
S}\neg\varphi_{k}\wedge\bigwedge_{k\in[m]\setminus S}\varphi_{k}\bigg{)}$
$\displaystyle=\sum_{S^{\prime}\supseteq S}(-1)^{|S^{\prime}\setminus
S|}\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}\bigg{(}\bigwedge_{k\in
S^{\prime}}\neg\varphi_{k}\bigg{)}$
$\displaystyle=\sum_{I\in\mathcal{I}(\Gamma):I\supseteq S}(-1)^{|I\setminus
S|}\prod_{k\in I}\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\neg\varphi_{k})$
$\displaystyle=\sum_{I\in\mathcal{I}(\Gamma):I\supseteq S}(-1)^{|I\setminus
S|}p_{I}=q_{S}.$

(When $S\notin\mathcal{I}(\Gamma)$, the above equality continues to hold, as
both side are zero.) Thus, the probability that clause $\varphi_{k}$ is false,
and all others true, is precisely $q_{\\{k\\}}$. The result now follows from
Theorem 5. ∎

Thanks to Lemma 2, the above results are completely robust against changes in
the implementation of algorithm PRS. Thus, the next scope to be resampled can
be selected by arbitrary means: the choice can be made on the current values
of variables, the past execution of the algorithm, or even externalities such
as random bits or the system clock. It is also valid to resample several
blocks simultaneously, in case several clauses are violated. If one is
interested in the expected number of individual variables resampled, this can
also be accessed though

(5) $\mathop{\mathbb{{}E}}\nolimits(\text{\\#variables
resampled})=\sum_{k=1}^{m}\frac{q_{\\{k\\}}a_{k}}{q_{\emptyset}},$

where $a_{k}$ is the arity of $\varphi_{k}$, for $1\leq k\leq m$.

##  3\. Example applications

One application, to sink-free orientations, will be done in detail to
illustrate the methods, and the other applications merely sketched. For ease
of presentation, all examples will be unweighted, i.e, the probability
distributions $\mathcal{D}_{i}$ are all uniform, as is the output
distribution. Incorporating weights does not require any conceptual changes.

###  3.1. Sink-free orientations of a graph

This approach to sampling sink-free orientations of a graph was introduced by
Cohn, Pemantle and Propp [4], and placed within the general framework of PRS
by Guo, Jerrum and Liu [14, §4.1].

Suppose $G=(V,E)$ is a graph with vertex set $\\{v_{1},v_{2},\ldots,v_{m}\\}$
and edge set $E=\\{e_{1},e_{2},\ldots,e_{n}\\}$.111The roles of $n$ and $m$
are reversed relative to the usual convention in graph theory, but this is
necessary to preserve consistency with the previous section. We wish to
sample, uniformly at random, an orientation of the edges of $G$ that has no
sinks, where a _sink_ is a vertex $v_{i}$ at which all incident edges are
oriented towards $v_{i}$. We assume that $G$ has at least one such sink-free
orientation. It is convenient to choose a reference orientation for the edges
of $G$ that is sink-free; denote by $\vec{G}$ the directed graph obtained from
$G$ by giving the edges of $G$ this reference orientation.

$v_{1}$$v_{2}$$v_{3}$$v_{4}$ $e_{1}$ $e_{2}$ $e_{3}$ $e_{4}$ $e_{5}$ Figure 3.
A sample graph $\vec{G}$ incorporating a reference orientation

To fit the pattern of PRS, we introduce Boolean variables
$X_{1},X_{2},\ldots,X_{n}$ and associate variable $X_{i}$ to edge $e_{i}$, for
$1\leq i\leq n$. These variables will be used to encode orientations of the
edges of $G$. The variable $X_{i}$ is to be interpreted as follows: if
$X_{i}=0$ then the edge $e_{i}$ is oriented against the reference orientation
(of $e_{i}$ in $\vec{G}$) and if $X_{i}=1$ then $e_{i}$ is oriented with the
reference orientation. Next, introduce clauses $\varphi_{k}$ to encode the
event that vertex $v_{k}$ is not a sink. So the scope of $\Phi_{k}$ is the set
${\sf Scp}(\Phi_{k})=\\{i:\text{$e_{i}$ is incident at $v_{k}$}\\}$, and the
clause $\varphi_{k}$ asserts that at least one edge incident at vertex $v_{k}$
is oriented away from $v_{k}$. By way of example, consider the graph $G$ in
Figure 3, which has been assigned a reference orientation to give a sink-free
directed graph $\vec{G}$. The condition that vertex $v_{2}$, for example, is
not a sink is asserted by the clause $\varphi_{2}(X_{1},X_{3},X_{4})=\neg
X_{1}\vee X_{3}\vee\neg X_{4}$. Then $\Phi$ is the formula

$\Phi(\mathbf{X})=(X_{1}\vee X_{2})\wedge(\neg X_{1}\vee X_{3}\vee\neg
X_{4})\wedge(\neg X_{2}\vee\neg X_{3}\vee X_{5})\wedge(X_{4}\vee\neg X_{5})$

that we encountered already in the previous section.

We observed earlier that $\Phi$ is an extremal instance. This is true in
general for sink-free orientations. If we have indices $1\leq k<\ell\leq m$
such that ${\sf Scp}(\varphi_{k})\cap{\sf Scp}(\varphi_{\ell})\not=\emptyset$
then necessarily vertices $v_{k}$ and $v_{\ell}$ are adjacent. But then it is
impossible for $v_{k}$ and $v_{\ell}$ to both be sinks, and hence
$\varphi_{k}\vee\varphi_{\ell}$ must hold. So Theorem 1 immediately assures us
that PRS will produce a uniform random sink-free orientation with probability
1. But is the expected running time polynomial in $n$ and $m$? In order to
apply Corollary 6 we need to bound the ratio $q_{\\{k\\}}/q_{\emptyset}$.
Although we don’t have a handle on $q_{\\{k\\}}$ and $q_{\emptyset}$ — and, in
a sense, $q_{\emptyset}$ is a quantity we would like to compute — we can bound
the ratio by defining an appropriate mapping from orientations with exactly
one sink, and those with none.

Introduce a function $f$ from $\\{1,\ldots,m\\}$ to itself that is consistent
with the reference orientation, that is to say, $(v_{k},v_{f(k)})$ is a
(directed) edge in $\vec{G}$ for all $1\leq k\leq m$. This is possible because
the reference orientation is sink-free. To each orientation of $G$ that has a
single sink at $v_{k}$ we associate a sink-free orientation as follows. Let
$e_{i}$ be the edge $(v_{k},v_{f(k)})$. Reverse the orientation of $e_{i}$,
i.e., set $X_{i}$, which was previously 0, to 1. Vertex $v_{k}$ is no longer a
sink, but $v_{f(k)}$ may have become one. If $v_{f(k)}$ is not a sink than
halt. Otherwise reverse the orientation of the edge $(v_{f(k)},v_{f^{2}(k)})$,
and continue. This process must terminate. For suppose not. Let $t$ be the
first instant at which we revisit a vertex, i.e., such that
$f^{t}(k)=f^{s}(k)$ for some $0\leq s<t$. The edge
$(v_{f^{s}(k)},v_{f^{s+1}(k)})$ is directed away from $v_{f^{s}(k)}$, and
hence vertex $v_{f^{s}(k)}$ is not a sink, a contradiction. (It is important
to note that we leave and revisit vertex $v_{s}$ via different edges.)

The edges that were flipped in the above construction form a path
$v_{k}=v_{f^{0}(k)},\allowbreak v_{f^{1}(k)},\ldots,\allowbreak
v_{f^{\ell}(k)}$. We may undo the construction provided we know $f^{0}(k)=k$
and $f^{\ell}(k)$. It follows that the number of orientations with a single
sink exceeds the number of sink free orientations by a factor at most $m^{2}$.
So by Corollary 6 the expected number of iterations in a run of PRS — in this
case the number of sinks that are ‘popped’ — is bounded above by
$m^{2}=|V(G)|^{2}$. We may also bound the number times the orientations of
individual edges are flipped. Fix a vertex $v_{k}$. We saw above how to repair
an orientation with a single sink at $v_{k}$. To undo this repair, we just
need to specify the index $f^{\ell}(k)$. Thus the number of orientations with
a single sink at $v_{k}$ exceeds the number of sink-free orientations by a
factor $m$. Referring to (5), we have $q_{\\{k\\}}/q_{\emptyset}\leq m$ and
$a_{k}=\deg(v_{k})$, the degree of vertex $v_{k}$. Thus the expected number of
edge orientation reversals is
$\sum_{k=1}^{m}q_{\\{k\\}}a_{k}/q_{\emptyset}\leq\sum_{k=1}^{m}m\deg(v_{k})\leq
2mn$. So the expected number of orientation reversals is a most
$2|V(G)|\,|E(G)|$. All this is in agreement with [4].

Surprisingly, the upper bound on edge-reversals can be tightened further to
$|E(G)|+|V(G)|^{2}$: see Guo and He [10]. Note that, the runtime analysis
critically used the assumption that coin tosses are unbiased, so that either
orientation of an edge is equally likely. (A simple counterexample shows that
this assumption is necessary.) In contrast, correctness of the algorithm
extends to asymmetric orientation probabilities.

###  3.2. Spanning trees of a graph

The Cycle-popping algorithm is an approach to uniformly sampling spanning
trees in a graph, introduced by Propp and Wilson [27]. Suppose $G$ is a graph
with vertex set $V=\\{v_{0},v_{1},\ldots,v_{n}\\}$ and edge set $E$. Instead
of sampling spanning trees in $G$ we will instead sample spanning
(in-)arborescences222That is, directed spanning trees with edges directed
towards a root vertex. rooted at $v_{0}$, which is of course equivalent.

For each $1\leq i\leq n$, define $D_{i}=\\{j:\\{v_{i},v_{j}\\}\in E\\}$, and
make $D_{i}$ into a probability space $\mathcal{D}_{i}$ by equipping it with
the uniform distribution. Introduce random variables $X_{1},\ldots,X_{n}$
distributed as $\mathcal{D}_{1},\ldots,\mathcal{D}_{n}$. These variables
indicate, for each $1\leq i\leq n$, a possible exit from vertex $v_{i}$. For
each simple (oriented) cycle
$C=(v_{i_{0}},v_{i_{1}},\ldots,v_{i_{\ell-1}},v_{i_{\ell}}=v_{i_{0}})$ define
the predicate $\varphi_{C}$ by

$\varphi_{C}=\neg(X_{i_{0}}=i_{1}\wedge X_{i_{1}}=i_{2}\wedge\cdots\wedge
X_{i_{\ell-1}}=i_{\ell}),$

and the formula $\Phi$ by $\Phi=\bigwedge_{C}\varphi_{C}$, where the
conjunction is over all oriented cycles in $G$. (In this context, ‘simple’ is
taken to mean ‘containing no repeated vertices’; thus we regard the 2-cycle
$(v_{i_{0}},v_{i_{1}},v_{i_{0}})$ as simple.) The intended interpretation of
the event $X_{i}=j$ is that vertex $v_{j}$ is the ancestor of vertex $v_{i}$
in the arborescence. The formula $\Phi(\mathbf{X})$ asserts that the ancestor
relation is consistent (has no cycles) and hence that $\mathbf{X}$ encodes a
spanning arborescence rooted at $v_{0}$.

Consider two clauses $\varphi_{C}$ and $\varphi_{C^{\prime}}$ corresponding to
distinct cycles $C$ and $C^{\prime}$. If ${\sf Scp}(\varphi_{C})\cap{\sf
Scp}(\varphi_{C^{\prime}})\not=\emptyset$ then $C$ and $C^{\prime}$ must have
a vertex in common. Select a vertex $v_{i}$ that is common to $C$ and
$C^{\prime}$ with the additional property that the successor to $v_{i}$ in
cycle $C$ is not equal to the successor to $v_{i}$ in cycle $C^{\prime}$. Let
$v_{j}$ be the successor to $v_{i}$ in $C$ and $v_{j^{\prime}}$ be the
successor in $C^{\prime}$. It is clear that $X_{i}=j$ and $X_{i}=j^{\prime}$
cannot both be true, and hence $\varphi_{C}$ and $\varphi_{C^{\prime}}$ cannot
both be false. Therefore $\Phi$ is extremal.

As, in the previous example, we need to estimate the ratio between
aborescences and ‘near-arborescences’ that contain a single cycle. (A near-
arborescence has two components: a spanning arborescence on some subset $S\ni
v_{0}$ of the vertices $V$ of $G$, rooted at $v_{0}$, and a unicyclic directed
subgraph spanning $V\setminus S$.) As before, by considering a suitable
mapping from near-arborescences to arboresecences, it can be shown that the
number of the former is at most $|V(G)|\,|E(G)|$ times the number of the
latter. Thus, by Corollary 6, the number of iterations made by PRS is at most
$|V(G)|\,|E(G)|$. A more refined analysis, due to Guo and He [10, Thm 15],
shows that the total number of _variable_ updates is bounded by essentially
the same expression.

###  3.3. Root-connected subgraphs

This ‘cluster-popping’ algorithm was proposed by Gorodezky and Pak [9], who
conjectured it to be efficient on a certain class of directed graphs. The
conjecture was resolved affirmatively by Guo and Jerrum [11].

Suppose $G=(V,A)$ is a directed graph with a distinguished root vertex $r\in
V$. A spanning subgraph $(V,S)$ of $G$ is said to be _root-connected_ if, for
every vertex $v\in V$, there is a directed path in $(V,S)$ from $v$ to $r$.
Our task is to sample, uniformly at random, a root-connected subgraph of $G$.
As usual, we restrict our attention to the unweighted version. However, as we
shall note later, the weighted version is of interest, owing to its connection
to a network reliability problem.

A subgraph $(V,S)$ may be encoded by variables $\mathbf{X}=(X_{e}:e\in A)$
taking values in $\\{0,1\\}$. The interpretation of $X_{e}=1$ is that $e\in
S$. For an arc $e\in A$, denote by $e^{-}$ and $e^{+}$ the start and end
vertex of $e$. A _cluster_ in $(V,S)$ is a set $\emptyset\subset C\subseteq
V\setminus\\{r\\}$ of vertices with the property that no edge $e\in S$ exists
with $e^{-}\in C$ and $e^{+}\in V\setminus C$. The property ‘$C$ is a cluster’
can be expressed formally by the predicate
$\psi_{C}=\bigwedge_{e\in\mathrm{cut}(C)}\neg X_{e}$, where
$\mathrm{cut}(C)=\\{e\in A:e^{-}\in C\text{ and }e^{+}\in V\setminus C\\}$. If
the subgraph $(V,S)$ has a cluster $C$ then it is clear that no vertex in $C$
can reach $r$, via a directed path in $(V,S)$, and hence $(V,S)$ is not root-
connected. The converse is also true: Suppose $(V,S)$ is not root-connected,
and let $v$ be some vertex from which the root $r$ is not reachable. Let $C$
be the set of all vertices reachable from $v$. Then $C$ is a cluster in
$(V,S)$.

This observation suggests that we should define

(6) $\Phi=\bigwedge_{\emptyset\subset C\subseteq
V\setminus\\{r\\}}\varphi_{C},$

where $\varphi_{C}=\neg\psi_{C}$. The formula $\Phi$ denies the existence of a
cluster in the subgraph encoded by $\mathbf{X}$, and hence correctly expresses
the property of being root-connected. The catch is that $\Phi$ is not in
general extremal. It is perfectly conceivable that two clusters $C,C^{\prime}$
exist that have nonempty intersection $C\cap C^{\prime}\not=\emptyset$. In
that case, we might have ${\sf Scp}(\varphi_{C})\cap{\sf
Scp}(\varphi_{C^{\prime}})\not=\emptyset$ and yet $\varphi_{C}$ and
$\varphi_{C^{\prime}}$ are both false. The solution is to make the predicates
$\varphi_{C}$ less demanding, while preserving the semantics of $\Phi$. We say
that the cluster $C$ is _minimal_ if it contains no cluster $C^{\prime}$ with
$C^{\prime}\subset C$. Then we define $\varphi_{C}$ to be true if $C$ is not a
minimal cluster. Formally,

$\varphi_{C}=\neg\Big{[}\psi_{C}\wedge\bigwedge_{\emptyset\subset
C^{\prime}\subset
C}\neg\psi_{C^{\prime}}\Big{]}=\neg\psi_{C}\vee\bigvee_{\emptyset\subset
C^{\prime}\subset C}\psi_{C^{\prime}}.$

Then define $\Phi$ as in (6). We claim that $\Phi$ still expresses the
condition that $\mathbf{X}$ encodes a root-connected subgraph $(V,S)$. If
$(V,S)$ is root-connected, then no cluster exists and hence $\varphi_{C}$ is
satisfied for all $\emptyset\subset C\subseteq V\setminus\\{r\\}$. Conversely,
suppose that $(V,S)$ is not root-connected. Then there is at least one
cluster, and hence at least one minimal cluster $C$. For this cluster,
$\varphi_{C}$ is contradicted, and hence $\Phi$ is false.

Although the meaning of $\Phi$ is unchanged, the formula is now extremal.
First note that, for all subsets $C$,

${\sf Scp}(\varphi_{C})=\\{X_{e}:e^{-}\in C\\}.$

So if $\varphi_{C}$ and $\varphi_{C^{\prime}}$ are any two distinct clauses
with ${\sf Scp}(\varphi_{C})\cap{\sf
Scp}(\varphi_{C^{\prime}})\not=\emptyset$, we must have $C\cap
C^{\prime}\not=\emptyset$. If $C$ and $C^{\prime}$ are both clusters then
$C\cap C^{\prime}$ must also be a cluster. Therefore, $C$ and $C^{\prime}$
cannot both be minimal clusters. It follows that at least one of $\varphi_{C}$
or $\varphi_{C^{\prime}}$ must hold. This deals with correctness of PRS in
this context.

Unfortunately, PRS does not have expected polynomial runtime on general
instances $G$, as can be appreciated by considering a counterexample presented
by Gorodezky and Pak [9]. However, those same authors conjectured that the
runtime is polynomial when the graph $G$ is ‘bidirected’, i.e., an edge exists
from vertex $u$ to $v$ in $G$ if and only if an edge exists from $v$ to $u$.
This special case is of interest, since root-connected subgraphs of a
bidirected graph $G$ correspond (via a constantly many-one relation) to
spanning connected subgraphs of the undirected version of $G$. Thus, cluster
popping provides a efficient approach to sampling connected spanning subgraphs
of a graph.

The conjecture of Gorodezky and Pak may be verified using Corollary 6. Again
the argument involves a mapping from subgraphs with exactly one minimal
cluster to root-connected subgraphs. The combinatorial details of this mapping
and its analysis, which are more involved in this case that the previous ones,
are given by Guo and Jerrum [11]. For a short while, PRS provided the only
known attack on sampling connected spanning subgraphs of a general undirected
graph, and its weighted version, undirected all-terminal reliability. However
the same problem (in a more general setting) has since been solved by Markov
chain simulation by Anari, Liu, Oveis Gharan and Vinzant [1].

###  3.4. Bases of bicircular matroids

Another application of PRS is to sampling bases of a bicircular matroid. The
algorithm was first presented in a slightly different guise by Kassel and
Kenyon [19]. Suppose $G=(V,E)$ is an undirected graph. The _bicircular
matroid_ associated with $G$ has $E$ as its ground set. The bases of the
matroid are all spanning subgraphs of $G$ in which every connected component
is unicyclic; equivalently, every connected component has the same number of
edges as it has vertices. The sampling algorithm may be derived methodically
using PRS. The application has similarities with the cycle-popping algorithm
described above in the context of sampling spanning trees.

As with cycle popping, variables are introduced that encode a function $g$
from $V$ to itself that respects the edges of $G$. (This is a slight deviation
from the spanning trees case, where the function was from $V\setminus\\{r\\}$
to $V$.) The spanning subgraph $(V,S)$ defined by $S=\\{\\{v,g(v)\\}:v\in
V\\}$ is very like a basis of the bicircular matroid, with two caveats. First,
we want to rule out cycles of length 2 — that is, situations in which
$g(g(v))=v$ for some $v\in V$ — as such functions $g$ do not correspond to
valid bases. Second, each basis with $c$ connected components corresponds to
$2^{c}$ distinct functions, as each cycle may be traced in either orientation.

To deal with these two objections, we specify a preferred orientation for
every cycle in $G$. Our formula $\Phi$ includes a clause $\varphi_{C}$, for
each potential cycle $C$ that either (a) has length two, or (b) is oriented in
the in the wrong sense. In each case, $\varphi_{C}$ asserts that $C$ does not
occur. It is easy to check that $\Phi$ is extremal. The expected number of
resamplings (either of clauses or individual variables) is $O(|V(G)|^{2})$.
Details are given by Guo and Jerrum [12].

###  3.5. Notes

The examples listed above are not the only known applications of PRS, but they
are the only non-trivial ones for which polynomial-time running time bounds
are known. At least, they are the one ones I am aware of.

One tempting potential application is to sampling strong orientations of an
undirected graph. An orientation of the edges of an undirected graph $G$ is
_strong_ if there is a directed path from every vertex of $G$ to every other.
If $G$ is connected, strong orientations coincide with ‘totally cyclic
orientations’. The number of total cyclic orientations of a graph $G$ is an
evaluation of the Tutte polynomial (at the point $(0,2)$). It is known that
counting totally cyclic (and hence strong orientations) is $\\#{\bf
P}$-complete [18]. However, the computational complexity of approximately
counting or uniformly sampling totally cyclic orientations is unknown.

The cluster-popping algorithm for root-connected subgraphs is easily adapted
to strong orientations. For a set $\emptyset\subset S\subset V=V(G)$ of
vertices of $G$, we say that $S$ is _cluster_ if all edges between $S$ and
$V\setminus S$ are directed into $S$. (The crucial difference with the root-
connected case is that there is no distinguished root vertex $r$ that is
excluded from all clusters.) We say that a cluster is _minimal_ if it is
minimal with respect to inclusion. As usual, define
$\Phi=\bigwedge_{\emptyset\subset S\subset V}\varphi_{S}$, where the formula
$\varphi_{S}$ expresses the condition that $S$ is not a minimal cluster. It
may be verified that $\Phi$ is extremal, and hence that PRS produces a uniform
random strong orientation (assuming that the $G$ has one, which happens
exactly when the graph $G$ is bridgeless). Unfortunately, the expected runtime
may be exponential, as can be appreciated by considering the ladder graph
$L_{n}$ on $2n$ vertices. (The ladder graph can be viewed as a $n\times 2$
rectangular piece of the square lattice, or as the cartesian product
$P_{n}\times P_{2}$ of a path on $n$ vertices and a path on $2$ vertices.) If
$v_{k}$ is one of the degree-2 corner vertices then the ratio
$q_{\\{k\\}}/q_{\emptyset}$ from Theorem 5 is exponential in $n$. (By
induction on $n$, the number of strong orientations is $2\times 3^{n-2}$,
whereas the number of orientations with a unique minimal cluster $\\{v_{k}\\}$
is at least $4^{n-1}$.) The fact that Theorem 5 gives an exact result and not
just an upper bound comes in useful here, as it enables us to deduce a _lower_
bound on the running time of PRS.

##  4\. Non-extremal instances

In an extremal instance, no two clauses that share variables can be
simultaneously false. We have seen that this leads to uniform outputs from
PRS. It transpires that we can get away with a little less than this. To this
end, introduce an explicit symmetric binary _dependency relation_
${\sim}\subseteq[m]^{2}$ on (indices) of clauses. If $k\sim\ell$ then the
scopes of clauses $\varphi_{k}$ and $\varphi_{\ell}$ intersect, i.e., ${\sf
Scp}(\varphi_{k})\cap{\sf Scp}(\varphi(\ell))\not=\emptyset$. This is as
before. However, the converse is not necessarily the case: there may be pairs
of clauses $\varphi_{k}$ and $\varphi_{\ell}$ that share variables for which
it is not the case that $k\sim\ell$. To get things to work we insist of the
following.

Axiom 1:

If $k\sim\ell$ then $\varphi_{k}\vee\varphi_{\ell}$ is a tautology.

Axiom 2:

Suppose $k\not\sim\ell$, $\neg\varphi_{k}$ and $\neg\varphi_{\ell}$. If we
resample the variables in ${\sf Scp}(\varphi_{\ell})$ then we are guaranteed
that $\neg\varphi_{k^{\prime}}$ for some clause $\varphi_{k^{\prime}}$ with
${\sf Scp}(\varphi_{k^{\prime}})\supseteq{\sf Scp}(\varphi_{k})$.

Note that an extremal instance satisfies the above axioms with $k^{\prime}=k$.
The additional flexibility allows PRS to be applied to a significantly wider
class of examples.

The algorithm is exactly as before except for one change. Although we have
considerable flexibility in the order in which to resample (the scopes of)
clauses, we no longer have complete freedom.

Algorithm 2 Partial Rejection Sampling with limited nondeterminism

$\mathrm{PRS}(\Phi,\mathcal{D})$ // $\Phi$ is a formula on variable set
$\mathbf{X}$

Sample $\mathbf{X}$ from the product distribution
$\mathcal{D}_{1}\times\cdots\times\mathcal{D}_{n}$

while $\neg\Phi(\mathbf{X})$ do

$N:=\\{\ell:\neg\varphi_{\ell}(X_{{\sf Scp}(\varphi_{\ell})})\\}$

Choose $k\in N$ deterministically, _based only on the set $N$ itself_

Resample all variables in ${\sf Scp}(\varphi_{k})$

end while

######  Theorem 7.

Suppose $\Phi$ is a satisfiable instance satisfying Axioms 1 and 2. Then
$\mathrm{PRS}(\Phi,\mathcal{D})$ terminates with probability 1. On
termination, $\mathbf{X}$ is a realisation of a random variable from the
distribution $\mathcal{D}_{\Phi}$.

###### Proof.

Termination with probability 1 can be argued exactly as in the proof of
Theorem 1.

For correctness, we set up the resampling table as in the proof of Theorem 1.
As before, fix a resampling table $T$, and run $\mathrm{PRS}$ on $T$ to obtain
a transcript. Since the algorithm has terminated, we know that the frontier
contains a satisfying assignment. Create a new resampling table $T^{\prime}$
by replacing the values in the frontier by some other satisfying assignment.
Now run the algorithm on the new resampling table $T^{\prime}$. We claim that
this second run correctly outputs the planted satisfying assignment.

If in both runs of the algorithm, the same clause $\varphi_{k}$ is selected in
every iteration then the output indeed will be correct. So assume that in some
iteration different clauses are selected in the two runs. Consider the first
iteration on which the set $N$ differs in the two runs, and suppose
$\varphi_{k}$ is true in one run and false in the other. As before, let
$S={\sf Scp}(\varphi_{k})$ and partition $S=S_{I}\cup S_{F}$, where variables
$S_{I}$ (respectively, $S_{F}$) take values from the interior (respectively,
frontier) of the resampling table. Note that $S_{I}\not=\emptyset$ (otherwise
$\varphi_{k}$ would be true in both runs) and $S_{F}\not=\emptyset$ (otherwise
$\varphi_{k}$ would have the same truth value in both runs).

There are two cases. Suppose first that $\varphi_{k}$ is false in the $T$-run
(and incidentally true in the $T^{\prime}$-run, thought this is not relevant
to the argument). Allow the $T$-run to continue. The algorithm cannot resample
${\sf Scp}(\varphi_{k})$ itself, since that action would take it past the
final frontier of the table. It cannot resample ${\sf Scp}(\varphi_{\ell})$
for any $\ell$ with $\ell\sim k$, by Axiom 1. It can resample ${\sf
Scp}(\varphi_{\ell})$ for some $\ell\not\sim k$. But, by Axiom 2, this action
would leave behind a clause $\varphi_{k^{\prime}}$ with $\varphi_{k^{\prime}}$
false and ${\sf Scp}(\varphi_{k^{\prime}})\supseteq{\sf Scp}(\varphi_{k})$. By
induction, at least one clause will always be false. Thus, we can never make
all clauses of $\Phi$ true, which contradicts the fact that the final frontier
contains a satisfying assignment.

The second and final case has $\varphi_{k}$ false in the $T^{\prime}$-run (and
incidentally true in the $T$-run). Up to this point, the two runs have made
exactly the same choices of scopes to resample. Now imagine that all the
resampling steps in the $T$-run are faithfully mirrored in the $T^{\prime}$
run. We have deviated from the deterministic choice rule of the algorithm, but
all resampling steps are legal, in the sense that we always resample scopes
${\sf Scp}(\varphi_{\ell})$ for which $\varphi_{\ell}$ is currently false. The
reason for this is exactly as in the proof of Theorem 1: briefly, that we
never resample variables in the frontier and the variables sampled from the
interior have the same values in both runs.

Since ${\sf Scp}(\varphi_{k})$ is never resampled again in the $T$-run, it is
also never resampled in the $T^{\prime}$-run. Also, we cannot resample a scope
${\sf Scp}(\varphi_{\ell})$ for any $\ell$ with $\ell\sim k$, by Axiom 1.
Finally, as in the first case, by resampling ${\sf Scp}(\varphi_{\ell})$ for
$\ell$ with $\ell\not\sim k$, we must leave behind a clause
$\varphi_{k^{\prime}}$ with $\varphi_{k^{\prime}}$ false and ${\sf
Scp}(\varphi_{k^{\prime}})\supseteq{\sf Scp}(\varphi_{k})$. Thus, we can never
make all clauses of $\Phi$ true, which contradicts the fact that the final
frontier contains a satisfying assignment.

Summarising, the final transcript remains unchanged if one satisfying
assignment is substituted for another in the final frontier. Thus, conditioned
on the transcript, each satisfying assignment
$\mathbf{X}=(b_{1},\ldots,b_{n})$ of $\Phi$ occurs with probability
proportional to
$\mathcal{D}_{1}(b_{1})\\*\mathcal{D}_{2}(b_{2})\cdots\mathcal{D}_{n}(b_{n})$.
So, at termination, $\mathbf{X}$ is distributed as $\mathcal{D}_{\Phi}$. ∎

###  4.1. Example: Independent sets

Suppose wish to sample independent sets in a graph $G$. Introduce variables
$\mathbf{X}=(X_{v}:v\in V(G))$ taking values in $\\{0,1\\}$ to encode
potential independent sets in $G$. The interpretation of $X_{v}=1$
(respectively $X_{v}=0$) is that vertex $v$ is in (respectively not in) the
independent set. In our product distribution we assume
$\mathop{\mathrm{Pr}}\nolimits(X_{v}=1)=\lambda/(1+\lambda)$ for all vertices
of $G$, for some positive ‘activity’ $\lambda$. (It is not essential that the
activity is constant over vertices, but it slightly simplifies the
exposition.) The natural formula expressing that $\mathbf{X}$ encodes an
independent set is $\Phi^{\prime}(\mathbf{X})=\bigwedge_{\\{u,v\\}\in
E(G)}\varphi_{\\{u,v\\}}$ where $\varphi_{\\{u,v\\}}=\neg(X_{u}\wedge X_{v})$.
However $\Phi^{\prime}$ is not extremal. Following the example provided by
cluster popping for root connected subgraphs, we try to re-express
$\Phi^{\prime}$ as a semantically equivalent extremal formula.

Let $S\subseteq V$ be a subset of at least two vertices that induces a
connected subgraph $G[S]$ of $G$. Denote by

$\partial S=\\{v:v\notin S\text{ and }\exists u\in S\text{ such that
}\\{u,v\\}\in E(G)\\}$

the boundary of $S$, containing all vertices outside of $S$ that are adjacent
to some vertex in $S$. We say that $S$ is a _cluster_ if $v\in S$ implies
$X_{v}=1$ and $v\in\partial S$ implies $X_{v}=0$. For each $S$ of the above
form, we introduce a clause $\varphi_{S}$ that asserts that $S$ is not a
cluster. Let $\Phi=\bigwedge_{S}\varphi_{S}$ where $S$ ranges over all vertex
subsets of size at least two that induce a connected subgraph. It is clear
that $\Phi(\mathbf{X})$ asserts that $\mathbf{X}$ encodes an independent set.

Unfortunately, a moment’s reflection reveals that $\Phi$ also is not extremal.
Denote by $\overline{S}$ the set $\overline{S}=S\cup\partial S$ and note that
${\sf Scp}(\varphi_{S})=\overline{S}$. It is possible to have clusters $S$ and
$S^{\prime}$ with $\overline{S}\cap\overline{S^{\prime}}\not=\emptyset$, in
which case ${\sf Scp}(\varphi_{S})\cap{\sf
Scp}(\varphi_{S^{\prime}})\not=\emptyset$, and yet
$\varphi_{S}\vee\varphi_{S^{\prime}}$ is false. For example, the path on five
vertices $\\{1,2,3,4,5\\}$ with $X_{3}=0$ and $X_{1}=X_{2}=X_{4}=X_{5}=1$ has
clusters $\\{1,2\\}$ and $\\{4,5\\}$, and
$\partial\\{1,2\\}\cap\partial\\{4,5\\}=\\{3\\}\not=\emptyset$.

However, suppose we define a dependency relation on clauses by $S\sim
S^{\prime}$ if and only if $S\cap\overline{S^{\prime}}\not=\emptyset$
(equivalently, if and only if $S^{\prime}\cap\overline{S}\not=\emptyset$).
Then it is straightforward to verify Axioms 1 and 2. If $S\sim S^{\prime}$
then $\varphi_{S}$ and $\varphi_{S^{\prime}}$ cannot be simultaneously false,
as that would entail setting each variable in $(S\cap\partial
S^{\prime})\cup(S^{\prime}\cap\partial S)\not=\emptyset$ both to 0 and to 1.
Also, if $S\not\sim S^{\prime}$ then resampling variables in
$\overline{S^{\prime}}=S^{\prime}\cup\partial S^{\prime}$ will leave the
variables in $S$ untouched. Either $\varphi_{S}$ will remain false, or some
clause $\varphi_{\widehat{S}}$ with $\widehat{S}\supset S$ will be made false.

Specialising the generic PRS algorithm to this example, we obtain the
following algorithm for sampling independent sets, which was first described
by Guo, Jerrum and Liu [14].

$\mathrm{PRSforIS}(G,\lambda)$ // $G$ is an undirected graph, and $\lambda$ a
positive real number

Sample $\mathbf{X}$ from the product distribution
$\mathrm{Bernoulli}(\lambda/(1+\lambda))^{n}$

while $\mathbf{X}$ does not encode an independent set do

Choose a cluster $S$ using a valid rule

Resample all variables $\\{X_{v}:v\in\overline{S}\\}$

end while

Algorithm 3 Partial Rejection Sampling for independent sets

######  Lemma 8.

$\mathrm{PRSforIS}(G,\lambda)$ terminates with probability 1. On termination,
$\mathbf{X}$ is a realisation of a random variable from the Gibbs distribution
for independent sets in $G$ with activity $\lambda$.

###### Proof.

Follows immediately from Theorem 7. ∎

###  4.2. Runtime Analysis

Sampling independent sets is in general an NP-hard problem [23, Thm 4], so we
need to make some assumption about the graph $G$ and activity $\lambda$. Our
goal is to find $\lambda_{\Delta}>0$ such that PRS terminates rapidly, for all
$\lambda<\lambda_{\Delta}$ and all graphs $G$ of maximum degree $\Delta$. We
no longer benefit from the tight analysis presented in Section 2 for extremal
instances. Let’s see how far we get using a simple analysis of how clusters
contract over time, on average. For convenience, identify the vertices of $G$
with the set $[n]$.

It is natural to measure the distance of a assignment $\mathbf{X}$ from
perfection by the integer
$D(\mathbf{X})=\sum_{S:\neg\varphi_{S}(\mathbf{X})}|S|$, in other words, the
number of variables set to 1 that find themselves in clusters. Suppose a
cluster $S$ with scope $\overline{S}$ is resampled. Let $c=|S|$ and
$b=|\partial S|$. Consider the profit/loss account for $D(\mathbf{X})$. Set
$z=\lambda/(1+\lambda)$. On the profit side, we lose the $c$ variables
$\mathbf{X}_{S}$, which are currently set to 1. After resampling, for any
given $i\in S$, the probability that variable $X_{i}$ contributes to
$D(\mathbf{X})$ is at most $\Delta z^{2}$. (This is a crude upper bound based
on the fact that $X_{i}$ and at least one of its neighbours must be assigned
1.) For a variable $X_{i}$ with $i\in\partial S$, the situation is worse.
Vertex $i$ may be adjacent to as many as $\Delta-1$ vertices $j$ lying outside
$\overline{S}$, each with $X_{j}=1$ and each surrounded by 0s. In this case,
with probability $z$, $D(\mathbf{X})$ increases by $\Delta$. (This is clearly
the worst case.) So the change in $D(\mathbf{X})$ is bounded above by $\Delta
cz^{2}+\Delta bz-c$. For fixed $c$, this expression is maximised when $S$
induces a tree in $G$, in which case $b=(\Delta-2)c+2$. The ratio $b/c$ is
maximised when that tree has 2 vertices. To be sure that $D(\mathbf{X})$
decreases in expectation we therefore require $\Delta
z^{2}+\Delta(\Delta-1)z-1\leq 0$. The critical value of $z$ at which equality
is achieved is given by

$z_{c}=\frac{1}{2}(\Delta-1)\left(\sqrt{1+\frac{4}{\Delta(\Delta-1)^{2}}}\,-1\right)$

which is asymptotic to $\Delta^{-2}$. Since $\lambda=z/(1-z)$, the critical
value $\lambda_{c}$ is also asymptotic to $\Delta^{-2}$. This is not the right
answer, as it is known that phase transition for hard-core model, both
computational and physical, is asymptotically $e/\Delta$. Although there is
scope for slightly tightening the above analysis, there is no hope of escaping
the $O(\Delta^{-2})$ asymptotics. The problem seems to be that we have to
assume the worst case environment for the cluster being resampled, a worst
case that is much worse than we would expect on average.

We start our search for an improved runtime bound with the analysis from
Section 2. One problem extending this analysis to the non-extremal situation
is that the proof of Theorem 5 fails. The reason for this is that the
interpretation of $q_{\emptyset}$ as the probability
$\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\Phi)$ that $\Phi$ is satisfied
is no longer valid. It transpires that this problem can be avoided by using a
different line of proof. Kolipaka and Szegedy [21, Thm 4] show that the number
of resamplings of the scope of $\varphi_{k}$ is bounded above by
$q_{\\{k\\}}/q_{\emptyset}$, provided the point $(p_{1},\ldots,p_{m})$ lies
within a certain region. This region was identified by Shearer [29] as the
theoretical limit of validity of the Lovász Local Lemma, even in the non-
algorithmic setting. Although elegant, it is difficult to use this result
directly: testing membership in the Shearer region in specific examples is
challenging, as is computing $q_{\\{k\\}}$ and $q_{\emptyset}$, which no
longer have simple combinatorial interpretations. Fortunately, there are
several weaker conditions that can be feasibly tested.

To progress, we need to assume the following condition.

Axiom 3:

The following is true relative to any resampling table, for any $k,\ell\in[m]$
with $k\not\sim\ell$. Suppose it is possible to resample ${\sf
Scp}(\varphi_{k})$ and then immediately resample ${\sf Scp}(\varphi_{\ell})$.
(Thus, $\neg\varphi_{k}$ initially, and after resampling ${\sf
Scp}(\varphi_{k})$ we have $\neg\varphi_{\ell}$.) Then it is possible to
resample ${\sf Scp}(\varphi_{\ell})$ and immediately resample ${\sf
Scp}(\varphi_{k})$.

It is easily checked Axiom 3 is satisfied by the independent sets example.
Consider two clauses $\varphi_{S}$ and $\varphi_{S^{\prime}}$ with $S\not\sim
S^{\prime}$. Suppose that the resampling block for ${\sf
Scp}(\varphi_{S^{\prime}})$ lies immediately above that for ${\sf
Scp}(\varphi_{S})$ in the resampling table. We observe that the order of the
resamplings can be reversed (leading to a locally different transcript). The
case ${\sf Scp}(\varphi_{S})\cap{\sf
Scp}(\varphi_{S^{\prime}})=\overline{S}\cap\overline{S^{\prime}}=\emptyset$ is
uninteresting. So consider a variable $X_{i}$ with
$i\in\overline{S}\cap\overline{S^{\prime}}$. Necessarily, $i\in\partial
S\cap\partial S^{\prime}$. It follows that $X_{i}$ takes the value 0 before
${\sf Scp}(\varphi_{S})$ is resampled (since $\neg\varphi_{S}$) and retains
that value after (since $\neg\varphi_{S^{\prime}}$). Therefore the two scopes
could as well have been resampled in the opposite order.

Axioms 1 and 3 imply that $([m],\sim)$ is a ‘lopsidependency graph’ [25, §6].
Axiom 3 is sometimes portrayed as as a positive dependency condition, but in
the resampling table view of the world it seems more natural to interpret it
as a commutativity condition. We say that a clause is _atomic_ if it is
falsified by exactly one assignment. A simpler but stronger alternative to
Axiom 3 is

Axiom $3^{\prime}$:

$\varphi_{k}$ is atomic, for every $k\in[m]$.

Axiom $3^{\prime}$ is even easier to check than Axiom 3. Again, Axioms 1 and
$3^{\prime}$ imply that $([m],\sim)$ is a lopsidependency graph.

Moser and Tardos [25, Thm 6.1] prove the following runtime bound.

######  Theorem 9.

Suppose that $\Phi$ is an instance with a dependency graph satisfying Axioms 1
and 3 (or Axioms 1 and $\textit{3\,}^{\prime}$). Suppose also that there
exists a sequence of reals $(x_{k}\in(0,1):k\in[m])$ such that, for all
$k\in[m]$,

$\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D}}(\neg\varphi_{k})\leq
x_{k}\prod_{\ell\in[m]:\ell\sim k}(1-x_{\ell}).$

Then, in expectation, Algorthm 2 resamples ${\sf Scp}(\varphi_{k})$ at most
$x_{k}/(1-x_{k})$ times before halting.

######  Lemma 10.

Suppose $G$ is a graph with $n$ vertices and maximum degree $\Delta$. There
exists $\lambda_{\Delta}>0$ such that the expected runtime of
$\mathrm{PRSforIS}(G,\lambda)$ is $O(n)$ whenever
$\lambda\leq\lambda_{\Delta}$. Asymptotically,
$\lambda_{\Delta}=\Theta(\Delta^{-1})$.

###### Proof.

Identifying vertices of $G$ with $[n]$, let

$\mathcal{C}=\big{\\{}S\subseteq[n]:|S|\geq 2\text{ and }G[S]\text{ is
connected}\big{\\}}$

be the set of all subsets of $V(G)$ that induce connected subgraphs of $G$
with at least two vertices. We need to find quantities
$\\{x_{S}:S\in\mathcal{C}\\}$ satisfying

$\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\lambda}(\neg\varphi_{S})\leq
x_{S}\prod_{S^{\prime}\in\mathcal{C}:S^{\prime}\sim S}(1-x_{S^{\prime}}).$

Note that we have included the activity $\lambda$ explicitly in the notation
here, as we are about to introduce a second artificial activity $\mu$. The
idea is to define the required quantities $x_{S}$ by
$x_{S}=\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\mu}(\neg\varphi_{S})$ for
some sufficiently small $\mu$ with $\mu>\lambda$. Suppose $S\in\mathcal{C}$
and let $c=|S|$ and $b=|\partial S|$. Also let $p=\mu/(1+\mu)$. We start with
a preliminary calculation. For any $i\in[n]$,

$\displaystyle\sum_{S^{\prime}\in\mathcal{C}:i\in S^{\prime}}x_{S^{\prime}}$
$\displaystyle=\sum_{S^{\prime}\in\mathcal{C}:i\in
S^{\prime}}\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\mu}(\neg\varphi_{S^{\prime}})$
(7)
$\displaystyle=\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\mu}\left(\bigvee_{S^{\prime}\in\mathcal{C}:i\in
S^{\prime}}\neg\varphi_{S^{\prime}}\right)$
$\displaystyle=\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\mu}(\text{$i$ is
contained in some cluster})$ (8) $\displaystyle\leq\Delta p^{2}.$

Equality (7) follows from disjointness of the events
$\neg\varphi_{S^{\prime}}$ over all $S^{\prime}\in\mathcal{C}$ with $i\in
S^{\prime}$ (Axiom 1). Inequality (8) is a simple upper bound on the
probability that vertex $i$ finds itself in a cluster, which we have already
seen. Then

$\displaystyle x_{S}\prod_{S^{\prime}\sim S}(1-x_{S^{\prime}})$
$\displaystyle\geq x_{S}\prod_{i\in S\cup\partial
S}\prod_{S^{\prime}\in\mathcal{C}:i\in S^{\prime}}(1-x_{S^{\prime}})$ (by over
counting) $\displaystyle\geq x_{S}\prod_{i\in S\cup\partial
S}\left(1-\sum_{S^{\prime}\in\mathcal{C}:i\in
S^{\prime}}x_{S^{\prime}}\right)$ $\displaystyle\geq p^{c}(1-p)^{b}(1-\Delta
p^{2})^{b+c}$ (by (8)) $\displaystyle\geq p^{c}(1-p)^{(\Delta-1)c}(1-\Delta
p^{2})^{\Delta c}$ (since $b\leq(\Delta-1)c)$) $\displaystyle=p_{\Delta}^{c},$
$\displaystyle p_{\Delta}$ $\displaystyle=p(1-p)^{\Delta-1}(1-\Delta
p^{2})^{\Delta}.$

Choose $p\in[0,1]$ to maximise $p_{\Delta}$, and note that both $p$ and
$p_{\Delta}$ are both $\Theta(\Delta^{-1})$. Let
$\lambda_{\Delta}=p_{\Delta}/(1-p_{\Delta})$. Then observe that for all
$\lambda\leq\lambda_{\Delta}$,

$\mathop{\mathrm{Pr}}\nolimits_{\mathcal{D},\lambda}(\neg\varphi_{S})\leq\bigg{(}\frac{\lambda_{\Delta}}{1+\lambda_{\Delta}}\bigg{)}^{c}=p_{\Delta}^{c}\leq
x_{S}\prod_{S^{\prime}\sim S}(1-x_{S^{\prime}})$

as required.

Having shown that the premise of Theorem 9 holds, we can read off an upper
bound on the expected number of resamplings. Repeating an earlier trick,

$\displaystyle\sum_{S\in\mathcal{C}}x_{S}$ $\displaystyle\leq\sum_{i\in
V(G)}\sum_{S^{\prime}\in\mathcal{C}:i\in S^{\prime}}x_{S^{\prime}}$
$\displaystyle\leq n\Delta p^{2}$ (by (8)) $\displaystyle=\Theta(n/\Delta).$

Since $x_{S}$ is bounded away from 1, the total number of resampled scopes is
also $\Theta(n/\Delta)$. ∎

######  Remark.

When $\Delta=3$, we find numerically that $\lambda_{3}=0.08797+$. There are
several steps of the calculation where slack is present, and could be reduced
at the expense of additional combinatorial complexity. With some extra
calculation, one could show that the expected size of a cluster containing a
given vertex $i$ is constant, and hence the expected number of _variable_
resamplings is linear in $n=|V(G)|$.

Extensive research on the algorithmic Lovaász Local Lemma has brought to light
a number of alternatives to Theorem 7. Examples which may be useful in
analysing PRS algorithms have been given by Bissacot, Fernández, Procacci, and
Scoppola [2], Harris [15], Harvey and Vondrak [16], Iliopoulos [17],
Kolmogorov [22], and Pegden [26]. A comprehensive treatment of the circle of
ideas surrounding the independent set polynomial $P_{\Gamma}$ and the Lovász
Local Lemma has been given by Scott and Sokal [28].

##  5\. Generalisations

In this article, we have restricted attention to the simplest version of PRS
based directly on the Moser Tardos algorithmic LLL. Specifically, we resample
the variables of just one clause at each step. This involved recasting the
‘obvious’ encoding of a problem as a CNF formula in a form suitable for
application of the method. For example, in the case of independent sets, we
replaced the natural two-variable clauses by larger clauses based on clusters.
Alternatively, it is possible to stick with the ‘natural’ formula at the
expense of complicating the resampling algorithm. This was the approach
originally taken by Guo, Jerrum and Liu [14]. Feng, Vishnoi and Yin [7]
extended the range of this approach by incorporating ‘soft’ constraints,
thereby allowing a wider range of spin systems to be addressed. Another
possible extension is to continuous state spaces, with Guo and Jerrum [13]
treating the hard-disks model, and Moka and Kroese [24] more general point
processes. Feng, Guo and Yin [6] show how to achieve perfect sampling when
strong spatial mixing holds. This last work is quite far from PRS, but still
relies only growing a sample by repeatedly ‘repairing’ parts of the current
configuration.

## Acknowledgements

The treatment of PRS presented here draws on many sources, in some cases
heavily. Particularly influential are the works of Moser and Tardos [25],
Knuth [20], Kolipaka and Szegedy [21] and Viennot [30]. I also learned a great
deal through collaboration with Heng Guo. Finally, in retrospect, it is
remarkable how many of the ideas behind PRS were already present in the work
of Propp and Wilson [27] on cycle-popping.

## References

  * [1] Nima Anari, Kuikui Liu, Shayan Oveis Gharan, and Cynthia Vinzant.  Log-concave polynomials II: High-dimensional walks and an FPRAS for counting bases of a matroid.  In STOC’19—Proceedings of the 51st Annual ACM SIGACT Symposium on Theory of Computing, pages 1–12. ACM, New York, 2019.
  * [2] Rodrigo Bissacot, Roberto Fernández, Aldo Procacci, and Benedetto Scoppola.  An improvement of the Lovász local lemma via cluster expansion.  Combin. Probab. Comput., 20(5):709–719, 2011.
  * [3] P. Cartier and D. Foata.  Problèmes combinatoires de commutation et réarrangements.  Lecture Notes in Mathematics, No. 85. Springer-Verlag, Berlin-New York, 1969.
  * [4] Henry Cohn, Robin Pemantle, and James G. Propp.  Generating a random sink-free orientation in quadratic time.  Electr. J. Comb., 9(1), 2002.
  * [5] Kimmo Eriksson.  Strong convergence and a game of numbers.  European J. Combin., 17(4):379–390, 1996.
  * [6] Weiming Feng, Heng Guo, and Yitong Yin.  Perfect sampling from spatial mixing.  CoRR, abs/1907.06033, 2019.
  * [7] Weiming Feng, Nisheeth K. Vishnoi, and Yitong Yin.  Dynamic sampling from graphical models.  SIAM J. Comput., 50(2):350–381, 2021.
  * [8] James Allen Fill and Mark Huber.  The randomness recycler: a new technique for perfect sampling.  In 41st Annual Symposium on Foundations of Computer Science (Redondo Beach, CA, 2000), pages 503–511. IEEE Comput. Soc. Press, Los Alamitos, CA, 2000.
  * [9] Igor Gorodezky and Igor Pak.  Generalized loop-erased random walks and approximate reachability.  Random Structures Algorithms, 44(2):201–223, 2014.
  * [10] Heng Guo and Kun He.  Tight bounds for popping algorithms.  Random Structures Algorithms, 57(2):371–392, 2020.
  * [11] Heng Guo and Mark Jerrum.  A polynomial-time approximation algorithm for all-terminal network reliability.  SIAM J. Comput., 48(3):964–978, 2019.
  * [12] Heng Guo and Mark Jerrum.  Approximately counting bases of bicircular matroids.  Combin. Probab. Comput., 30(1):124–135, 2021.
  * [13] Heng Guo and Mark Jerrum.  Perfect simulation of the hard disks model by partial rejection sampling.  Ann. Inst. Henri Poincaré D, 8(2):159–177, 2021.
  * [14] Heng Guo, Mark Jerrum, and Jingcheng Liu.  Uniform sampling through the Lovász local lemma.  J. ACM, 66(3):Art. 18, 31, 2019.
  * [15] David G. Harris.  Lopsidependency in the Moser-Tardos framework: beyond the lopsided Lovász local lemma.  ACM Trans. Algorithms, 13(1):Art. 17, 26, 2016.
  * [16] Nicholas J. A. Harvey and Jan Vondrák.  An algorithmic proof of the Lovász local lemma via resampling oracles.  SIAM J. Comput., 49(2):394–428, 2020.
  * [17] Fotis Iliopoulos.  Commutative algorithms approximate the LLL-distribution.  In Approximation, randomization, and combinatorial optimization. Algorithms and techniques, volume 116 of LIPIcs. Leibniz Int. Proc. Inform., pages Art. No. 44, 20. Schloss Dagstuhl. Leibniz-Zent. Inform., Wadern, 2018.
  * [18] F. Jaeger, D. L. Vertigan, and D. J. A. Welsh.  On the computational complexity of the Jones and Tutte polynomials.  Math. Proc. Cambridge Philos. Soc., 108(1):35–53, 1990.
  * [19] Adrien Kassel and Richard Kenyon.  Random curves on surfaces induced from the Laplacian determinant.  Ann. Probab., 45(2):932–964, 2017.
  * [20] Donald E. Knuth.  The Art of Computer Programming, volume 4, Fascicle 6.  Addison-Wesley Professional, 2015.
  * [21] Kashyap Kolipaka and Mario Szegedy.  Moser and Tardos meet Lovász.  In STOC’11—Proceedings of the 43rd ACM Symposium on Theory of Computing, pages 235–243. ACM, New York, 2011.
  * [22] Vladimir Kolmogorov.  Commutativity in the algorithmic Lovász local lemma.  SIAM J. Comput., 47(6):2029–2056, 2018.
  * [23] Michael Luby and Eric Vigoda.  Fast convergence of the Glauber dynamics for sampling independent sets.  Random Structures Algorithms, 15(3-4):229–241, 1999.
  * [24] Sarat B. Moka and Dirk P. Kroese.  Perfect sampling for Gibbs point processes using partial rejection sampling.  Bernoulli, 26(3):2082–2104, 2020.
  * [25] Robin A. Moser and Gábor Tardos.  A constructive proof of the general Lovász Local Lemma.  J. ACM, 57(2), 2010.
  * [26] Wesley Pegden.  An extension of the Moser-Tardos algorithmic local lemma.  SIAM J. Discrete Math., 28(2):911–917, 2014.
  * [27] James G. Propp and David B. Wilson.  How to get a perfectly random sample from a generic Markov chain and generate a random spanning tree of a directed graph.  J. Algorithms, 27(2):170–217, 1998.
  * [28] Alexander D. Scott and Alan D. Sokal.  The repulsive lattice gas, the independent-set polynomial, and the Lovász local lemma.  J. Stat. Phys., 118(5-6):1151–1261, 2005.
  * [29] James B. Shearer.  On a problem of Spencer.  Combinatorica, 5(3):241–245, 1985.
  * [30] Gérard Xavier Viennot.  Heaps of pieces. I. Basic definitions and combinatorial lemmas.  In Combinatoire énumérative (Montreal, Que., 1985/Quebec, Que., 1985), volume 1234 of Lecture Notes in Math., pages 321–350. Springer, Berlin, 1986.


¡html¿¡head¿ ¡meta http-equiv=”content-type” content=”text/html;
charset=ISO-8859-1”¿

¡title¿CERN-2014-009¡/title¿

¡/head¿

¡body¿

¡h1¿¡a href=”http://cas.web.cern.ch/cas/Norway-2013/Norway-advert.html”¿CAS -
CERN Accelerator School: Advanced Accelerator Physics Course (CAS
2013)¡/a¿¡/h1¿

¡h2¿Trondheim, Norway, 18-29 Aug 2013¡/h2¿

¡h2¿Proceedings - CERN Yellow Report

¡a href=”https://cds.cern.ch/record/1507631?ln=en”¿CERN-2014-009¡/a¿¡/h2¿

¡h3¿Editors: W. Herr¡/h3¿

This report presents the proceedings of the Course on Advanced Accelerator
Physics organized by the CERN Accelerator School. The course was held in
Trondheim, Norway from 18 to 29 August 2013, in collaboration with the
Norwegian University of Science and Technology. Its syllabus was based on
previous courses and in particular on the course held in Berlin 2003 whose
proceedings were published as CERN Yellow Report CERN- 2006-002. The field has
seen significant advances in recent years and some topics were presented in a
new way and other topics were added. The lectures were supplemented with
tutorials on key topics and 14 hours of hands on courses on Optics Design and
Corrections, RF Measurement Techniques and Beam Instrumentation and
Diagnostics. These courses are a key element of the Advanced Level Course.
¡h2¿Lectures¡/h2¿

¡p¿ LIST:arXiv:1601.04901¡br¿ LIST:arXiv:1601.04907¡br¿
LIST:arXiv:1601.04913¡br¿

¡p¿ Title: ¡a href=”https://cds.cern.ch/record/1982420?ln=en”¿RF Measurement
Concepts¡/a¿ ¡br¿ Author: F. Caspers and P. Kowina¡br¿ ¡/p¿

LIST:arXiv:1601.05411¡br¿ LIST:arXiv:1601.05201¡br¿ LIST:arXiv:1601.05207¡br¿
LIST:arXiv:1507.02213¡br¿ LIST:arXiv:1601.05210¡br¿ LIST:arXiv:1601.05214¡br¿
LIST:arXiv:1601.05216¡br¿ LIST:arXiv:1601.05227¡br¿ LIST:arXiv:1601.05230¡br¿
LIST:arXiv:1601.05235¡br¿ ¡p¿ Title: ¡a
href=”https://cds.cern.ch/record/1982431?ln=en”¿High-Field Accelerator
Magnets¡/a¿ ¡br¿ Author: G. de Rijk¡br¿ ¡/p¿ LIST:arXiv:1601.05259¡br¿
LIST:arXiv:1601.05258¡br¿ ¡/p¿ ¡/body¿¡/html¿