{'fraction_non_alphanumeric': 0.03149088736181655, 'fraction_numerical': 0.010457125784284434, 'mean_word_length': 4.494418910045962, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Exact sampling from non-attractive distributions using summary states

Propp and Wilson's method of coupling from the past allows one to efficiently generate exact samples from attractive statistical distributions (e.g., the ferromagnetic Ising model). This method may be generalized to non-attractive distributions by the use of summary states, as first described by Huber. Using this method, we present exact samples from a frustrated antiferromagnetic triangular Ising model and the antiferromagnetic q=3 Potts model. We discuss the advantages and limitations of the method of summary states for practical sampling, paying particular attention to the slowing down of the algorithm at low temperature. In particular, we show that such a slowing down can occur in the absence of a physical phase transition.


I. INTRODUCTION
In many statistical problems, physical and otherwise, it is useful to be able draw samples from a complex distribution. For example, in statistical physics one is interested in the Boltzmann distribution where E(σ) describes the energy of a system in configuration σ, β is the inverse temperature (we set k B = 1), and Z is a normalizing constant (the partition function). In general, E(σ) may be easy to evaluate for a particular configuration, but the number of possible configurations makes it impractical to draw directly from the distribution. Yet some efficient method of sampling is desirable, as this would allow one to calculate properties of the system that might not be easily computed by analytical means.
In traditional Monte Carlo sampling methods, such as the Metropolis-Hastings method [1] and Gibbs sampling [2] (also known as the heat bath algorithm), one constructs an ergodic Markov chain whose stationary distribution is the desired distribution. By starting in some state and evolving the chain for a sufficiently long time, one can approximate a sample from the desired distribution. Unfortunately, such a sample is exact only in the limit of infinite time. In practice, it is often difficult to determine how long to wait to achieve sufficiently good samples, and one inevitably either produces poor samples or wastes time by running the Markov chain for longer than necessary.
However, in 1996, Propp and Wilson demonstrated the possibility of exact sampling by the method of coupling from the past, allowing one to produce perfect samples in a finite number of steps [3]. In the most general case, their method requires the infeasible task of running a Markov chain for every possible initial state of the system. But for certain distributions, termed attractive (such as a ferromagnetic Ising model), Propp and Wilson showed that the task may be greatly simplified by tracking only extremal states, permitting the practical calculation of exact samples. This method was generalized to anti-attractive distributions by Häggström and Nelander [4]. More recently, Huber showed that one may instead track only a single state that summarizes one's knowledge of the system [5]. Because it does not require that the states be partially ordered, this last method is applicable to non-attractive distributions.
Using the summary state method, we have drawn exact samples from the antiferromagnetic triangular Ising model and from the Potts model. Fig. 1 shows one such sample. In §II, we describe the methods that make this possible. In §III, we briefly discuss the Ising and Potts models. We present results from the exact sampling of these models in §IV. Finally, we discuss the convergence properties of the summary state method, and we suggest practical generalizations. Propp and Wilson's method of coupling from the past [3] is based on the observation that, for a fixed choice of the random numbers used to propagate a Markov chain, its possible paths in state space may ultimately coalesce into a single trajectory. Once two initial states lead to the same state, they will remain in lock step thereafter.
Consider simulating a Markov chain from every possible initial state at some fixed time t = −T , with the goal of taking a sample at t = 0. If all the chains coalesce before t = 0, then this finite procedure yields the same results as a Monte Carlo simulation started at an infinite time in the past, so the result is an exact sample. If the chains fail to coalesce, one can simply double the starting time to −2T , reusing the random numbers for the interval [−t, 0] (i.e., treating the random numbers as a function of simulation time), and repeat until coalescence is achieved.
Having to follow every possible state would make this method exponentially intractable. But for problems that admit a partial ordering of the states and which are "attractive" -that is, which preserve the ordering under evolution of the Markov chain -the computation can be vastly simplified by tracking only the extremal states. An example of an attractive system is the ferromagnetic Ising model, in which it is energetically favorable for spins to align with each other.
Huber [5] and Harvey and Neal [6] have shown that the method of Propp and Wilson may be generalized using a single summary state instead of a pair of extremal states. This single state summarizes one's knowledge of the possible states of the system, allowing the state of some subsystems to be uncertain.
For example, suppose the system is a collection of variables σ i taking on the values {±1}. Conventional Gibbs updating sets where u is uniformly distributed on [0, 1] andσ i denotes the set of all variables but the ith. To implement summary states, we allow each variable to take on the additional value ? which indicates uncertainty. We then run a modified Markov chain on this system: σ i is updated according to Eq. (2) if the result is the same for any possible assignment of ±1 to the ?'s inσ i ; otherwise, σ i → ?. As in the Propp and Wilson method, we run the chain from successively longer times in the past with random numbers as a function of simulation time. When no variables remain in ? states, the algorithm has converged, and we may take a sample at t = 0. For the case of attractive distributions, this procedure is exactly equivalent to the Propp and Wilson scheme. The value ? denotes variables that differ between the maximal and minimal states, and removal of all ? states corresponds to coalescence of the bounding chains. However, using a single summary state, there is no requirement that the states be ordered in any way. Thus the summary state method can also be applied to nonattractive distributions -for example, the antiferromagnetic Ising model.
Although the samples returned by this method are exact, the algorithm does not necessarily converge after a reasonable amount of time. Huber has shown that for antiferromagnetic spin systems at sufficiently high temperature, the expected running time of the algorithm is polynomial in the number of spins [5]. However, for systems with a phase transition, the convergence time diverges as a power law at the critical temperature, a phenomenon known as critical slowing down [7].

Consider the Hamiltonian
where J mn is the coupling between spins m and n and H m is the value of an external magnetic field at the location of spin m. The appropriate Markov chain update rule is Eq. (2) with In the Ising model [8], J mn is taken to be zero unless spins m and n are adjacent, in which case it is some constant J. Cases of particular interest are the square lattice, in which each spin has four neighbors, and the triangular lattice, with six neighbors per spin. In general, the behavior of Ising systems can vary with their spin connectivity. For both kinds of lattices, we use periodic boundary conditions. Because we may simultaneously update the states of spins whose conditional distributions are independent, one iteration of the Markov chain consists of two sweeps for the square lattice and three for the triangular lattice.
In this paper, we use the normalization J = ±1. J = +1 corresponds to the ferromagnetic case, in which spins prefer to point in the same direction; J = −1 corresponds to the antiferromagnet. As mentioned previously, the ferromagnetic case is attractive. The antiferromagnet on a square lattice is a special case, because its properties are isomorphic to those of a square ferromagnet. However, for a triangular lattice, there is no such isomorphism. With six neighbors per spin, there is no way to minimize the energy locally at all sites: we say the system is frustrated.
It is well known that a two-dimensional ferromagnetic Ising model exhibits a phase transition [9]. Below a critical temperature β −1 c , there is spontaneous symmetry breaking, and the system develops a preferred spin orientation in the absence of any magnetic field. For a square lattice, β −1 c = 2.27. At this temperature, the relaxation time of the dynamic system diverges, a phenomenon known as critical slowing down [7]. Correspondingly, there is a divergence in the convergence time for some Markov chain Monte Carlo algorithms, such as coupling from the past, and exact samples cannot be generated for lower temperatures. Note that there is no phase transition in the case of a triangular antiferromagnet [10], so there cannot be a critical slowing down in the traditional sense.
To circumvent the problem of nonconvergence below the critical temperature, Propp and Wilson actually used a related system, the random cluster model, to generate Ising samples [3]. Unfortunately, this model has no obvious analog in the antiferromagnetic case.
The Potts model is a generalization of the Ising model wherein spins may take on q different values {0, 1, . . . , q − 1} [11,12]. Spins interact only with others of the same type. The Hamiltonian is Specifically, we consider the antiferromagnetic Potts model with q = 3 on a square lattice with zero magnetic field. This model has a critical point only at β −1 = 0, so there is no phase transition [13].

A. Ising model
By implementing the summary state method, we have produced exact samples from the Ising and Potts models. For example, Fig. 1 shows a sample from a triangular Ising antiferromagnet consisting of 120 2 = 14, 400 spins at β −1 = 4.9 with zero applied magnetic field.
We find that the number of iterations required for the algorithm to converge diverges at a threshold temperature. We have studied this divergence using a lattice of N = 63 2 = 3969 spins. Simulations using larger N (e.g., N = 99 2 ) suggest that the outcome is not significantly affected by choosing a larger grid size. Fig. 2 shows the divergence, to which we have fitted a power law of the form We find that the time diverges with an exponent b = 1.03 ± 0.01 at the threshold temperature β −1 t = 4.839 ± 0.005. This divergence is an important feature of the summary state method. It is qualitatively similar to critical slowing down, but note that no physical phase transition is involved. In divergent situations the augmented Markov chain has a metastable set of distributions with many ?'s, such that it is very unlikely for it to enter a state with no ?'s.
To draw an exact sample using the summary state method, the system must go from a completely uncertain state to a completely certain state. Thus, it must pass through a state with only a few scattered ?'s. For temperatures sufficientlty near the threshold, where we know that such a sparse configuration can be reached, we might expect that the limiting factor is the probability that an isolated ? can cause divergence.
Therefore, as a very rough estimate, we might suppose that the divergence occurs when the probability of a single ? turning one of its six neighbors into a ? rises above 1 6 . We expect that the neighbors of any given spin σ i should be (on average) half up and half down. Replacing one of these neighbors by a ?, we may assume the configuration (↑↑↑↓↓ ?) without loss of generality. Then the threshold temperature is determined by which has the solution β −1 = 4/ ln 2 ≈ 5.8.
To examine the validity of a threshold temperature analysis based on the persistence of single ?'s, we compiled statistics on the stability of an equilibrium system with a single ? added. Because we cannot create exact samples for much of the temperature range of interest, we generated approximate samples by simulating for fixed time (100 iterations) a random initial state. We then set one spin to ? and simulated the system forward. If any uncertainty remained after 500 iterations, we said the system diverged. Fig. 3 shows the fraction of divergent trials for various temperatures. As one would expect, this fraction goes to zero very near the threshold temperature. It is also interesting to consider how the algorithm behaves when a uniform nonzero magnetic field H is applied. Biasing the spins makes it easier for them to choose a particular orientation, so we would expect convergence to be easier. Fig. 4 shows the region of convergence in the (β −1 , H) plane.

B. Potts model
In addition, we have implemented exact sampling of the Potts model for arbitrary q. Fig. 5 shows an exact sample with q = 3 for a square antiferromagnetic lattice of 100 2 = 10, 000 spins. A naïve implementation of the summary state method would augment the possible spin values with a single ?. We refer to this method as algorithm A. However, it is possible to retain more information about uncertain spins: for each spin, we store a binary q-bit vector (b 1 , b 2 , . . . , b q ), b i ∈ {0, 1}. Bit b i is set to one if it is possible for the spin to take on the value i: thus the initial state of each spin is b = (1, 1, . . . , 1). In updating the state of the system, we set b i = 0 only when the spin cannot take on the value i for any allowed configuration of its neighbors. We refer to the latter method as algorithm B.
To demonstrate the advantage of retaining more information in the summary state, we have studied the convergence properties of both algorithms. This comparison is shown in Fig. 6, based on data for a square 64 2 = 4096 spin lattice. As in the Ising study, both algorithms lead to a power law divergence with an exponent of one (b A = 1.04 ± 0.03, b B = 0.99 ± 0.02). However, the threshold temperatures for the two algorithms are quite different: β −1 t,A = 2.293 ± 0.005, whereas β −1 t,B = 1.157 ± 0.004. As in the Ising example above, neither of the divergences corresponds to a physical phase transition.

V. CONCLUSIONS
We have demonstrated the usefulness of the summary state method for exact sampling from non-attractive distributions. In both the antiferromagnetic Ising and Potts models, the method works above a certain threshold temperature, with a power law divergence in the coalescence time at the threshold. Although similar to the phenomenon of critical slowing down, this divergence does not occur at a physical phase transition. As the Potts example shows, the location of the divergence is a feature of the specific implementation of the summary states, not of the underlying distribution. We have shown that retaining more information in the summary state will allow convergence at lower temperatures.
Based on this result, we may propose an improved algorithm for the triangular antiferromagnetic Ising model. At lower temperatures, the system should be increasingly ordered, and tracking this order might make it easier to gain incremental knowledge of the state of the system. One idea is to keep track of correlations between spins by grouping them into hexagonal clumps of seven, which can be used to tile the triangular lattice. Each tile has 2 7 = 128 possible states. In analogy to the Potts method presented earlier (in which a q-bit vector represents the uncertainty about a spin), representing each tile with a 128-bit vector would allow individually tracking the possible arrangements of those seven spins. Within a tile, the summary state can track anticorrelation, which we expect to arise at low temperatures. Each edge of a tile can be easily summarized in a 4-bit vector for comparison with its neighboring tiles. Each tile can then update its summary state by considering the possible states of the neighboring edges.
Of course, summary state sampling remains a valuable tool even if it cannot be done below some threshold. Where it does work, it is exact. The method also has the advantage that we need not even know if the distribution in question is attractive -it can be applied to any system.

VI. ACKNOWLEDGMENTS
We wish to thank Radford Neal and David Wilson for several helpful discussions. AMC and RBP acknowledge the support of the Caltech Cambridge Scholars Program. DJCM's group is supported by the Gatsby Charitable Foundation and by a Partnership Award from IBM Zürich.