Metropolis, simulated annealing and proposal distributions

Question

I'm trying to understand the physical meaning of the proposal step generating function in Metropolis algorithm.

In the original paper, and most derivations I found, it seems that it's not much important as long as it's symmetric (for detailed-balance) and allows a good exploration of the configuration space (e.g. $\sim 0.5$ acceptance ratio). Usually it's either a uniform or a normal distribution.

Now when studying Simulated Annealing, most papers assume a common use of a normal distribution with $\sigma = T$.

Given simulated annealing is deeply rooted in physical analogy, if we put the temperature there, there must be some physical sense we're trying to model. I mean, we could use a fixed uniform distribution with the proper width and still get good enough parameter space exploration and boltzmann distribution in the result, so why we use a gaussian and why the temperature dependence?

Also, when Cauchy distribution is introduced (e.g. by Szu and Hartley, 1987) they say the temperature used as a control parameter is kind of artificial, implying that the one in the Gaussian was not.

So what's the physical phenomenon with Gaussian, temperature dependent, steps we're modeling?

I know it all can be better explained with markov chains, priors, posteriors and likelihood, but that way you lose the physical analogy...

Answer 1

I suppose one could draw a physical analogy, noting that the proposal length scale should be proportional to the kinetic energy of the particle, but as @Xi'an I too would argue that this is simply a practical choice to maintain sensible acceptance rates as the temperature goes down.

To see this, assume we want to optimize a typical statistical problem, where the target function is a likelihood, which is asymptotically normal. Assume for simplicity that the optimum is at 0, and the width (sd) is $\sigma_t$. We optimize the log Likelihood, which is proportional to

$$ LL \propto - x^2 / \sigma_t^2 $$

x being the distance from the maximum at 0. In the standard simulated annealing, our acceptance probability would thus be proportional to

$$ p_a \propto e^{- \frac{x^2} {\sigma_t^2 \cdot T } } $$

Now, consider that the expected value of the random variable x depends on the proposal function

$$ <x^2> \sim proposal $$

As this is tagged self-study, I'll leave it to you to see how you should scale the proposal function to get constant acceptance rates (I would argue that this is a sensible goal).

Answer 2

Here's an analogy from physics - consider the distribution in $n$ dimensions to have many peaks or valleys throughout. In a MCMC sampling problem, the 'height' of the peaks would be how improbable the state was - we would want to go to the valleys. Now, let us imagine this same surface to be a the map of potential energy for a physical system. Then, since a particle or object would want to find a minimum for the potential, it will try to go 'downhill'.

Now, in a classical sense, the particle can never pass through or 'tunnel' through the potential barriers. These barriers are the peaks in our distribution. However, in quantum theory, these particles always have a nonzero probability of tunneling through these peaks, and the probability of this tunneling behavior is dependent on the size of the barrier and the amount of energy the particle has.

Essentially, if the particle has a lot of energy, it is moving fast. If it is moving fast, it can tunnel through the potential barrier. If a particle's aggregate phase of matter has a high temperature, each particle has a lot of kinetic energy, on average. Let's reconsider our analogy with simulated annealing. If we have particles with high temperature, we allow them to escape local minima by going through otherwise unfavorable portions of state space. They are, in this analogy, 'tunneling' through the regions of low probability.

Now, as to how this relates to the distributions, I'm not sure.

Answer 3

In the original paper, and most derivations I found, it seems that it's not much important as long as it's symmetric (for detailed-balance) and allows a good exploration of the configuration space (e.g. ∼0.5

acceptance ratio). Usually it's either a uniform or a normal distribution.

Now when studying Simulated Annealing, most papers assume a common use of a normal distribution with σ=T .From a physicist point of view simulated annealing seeks the most probable state of a system, described by the Boltzmann distribution: $$ w(p,q|T) \propto \exp\left( -\frac{H(p,q)}{T}\right) $$ where $H(p,q)$ is the Hamiltonian, also referred to as the energy of the system, whereas $T$ is temperature. Generalized momenta and coordinates, $p$ and $q$ are related via Hamiltonian evolution (heavily exploited, e.g., in the Hamiltonian Monte Carlo); these need not he continuous - e.g., they could by discrete up/down states of magnetic moments in ferromagnet or placement of atoms A and B in a binary alloy. MCMC is used to sample different variable configurations not connected via Hamiltonian evolution - physically it mostly means sampling states of different energies. This does not imply a gaussian proposal density and largely depends on the nature of the variables and the Hamiltonian.

The essential part is that at the temperature of interest, $T_0$ the system us likely to be found in its lowest energy state. Finding this state corresponds to maximizing the probability/likelihood $w(p,q|T_0)$.

In case of a multimodal likelihood we are not guaranteed to find the global maximum - physically it means that there may exist multiple global maxima (e g., corresponding to different directions of magnetization in a ferromagnet) or there may be some metastable/long-living "glass" states, the transitions from which to the global maximum are unlikely, since they require passing through highly improbable intermediate configurations.

As the probability of a state is proportional to $\exp\left( -\frac{H(p,q)}{T}\right)$, intermediate states are easier reached, if the temperature $T$ is higher. The solution is then to explore the energy states at a high temperature, and then use the result to redo the exploration at lower and lower temperature, till we reach the desired one.

Probabilistically, instead of exploring distribution $$p(x)\propto \exp\left(-l(x)\right),$$ we work with a rescaled negative log-probability $$p_\alpha(x)\propto \exp\left(-\alpha l(x)\right),\alpha \leq 1,$$ which is more easily amenable to exploration, e.g., via the Metropolis-Hastings algorithm.

One particular class of physical problems is particles undergoing diffusive/Langevin dynamic. These can be modeled as a random walk with a gaussian proposal density. They are referred to sometimes as drift-diffusion problems, but they are rather ubiquitous in science and may appear under many different names. The easy "physical" way of thinking of such systems is water on a rugged surface, which is being shaked. When the holes in the surface are deep, while the shaking is weak, the water remains where it is - the shaking explores only the surface near the minima. Shaking harder would allow the water to leak to other places. Annealing achieves it by rescaling the rugged landscape, allowing the water flow, and then gradually restoring the original size of the roughness.

Update

In the original paper, and most derivations I found, it seems that it's not much important as long as it's symmetric (for detailed-balance) and allows a good exploration of the configuration space (e.g. ∼0.5 acceptance ratio). Usually it's either a uniform or a normal distribution.

Now when studying Simulated Annealing, most papers assume a common use of a normal distribution with σ=T

Near the probability maximum (energy minimum) log-likelihood can be considered as parabolic, i.e., we have $$ p(x)\propto \exp\left(-\frac{kx^2}{T}\right), $$ which means that the characteristic extent of the distribution is $\sigma=\sqrt{T/k}$. This gives the idea for the suitable step of the Metropolis chain to get the acceptance rate around 0.5: if the proposal distribution is much wider, most steps will be rejected, while if it is much narrower, most of them will be accepted.

As it is seen from the above argument, we are not talking here literally about $\sigma=\sqrt{T}$, but rather about the proposal distribution width scaling with temperature as $\sqrt{T}$ - there is still a factor that depends on the curvature of the probability distribution near its peak.

---

References

1. Simulated annealing: Theory and applications (1987) by Laarhoven and Aarts is an early book taking simulated annealing from physics domain to general statistical applications.
2. Reaction rate: fifty years after Kramers (1990) by Hänggi et al. is a review article describing many physical, chemical and other systems, where diffusion-like dynamics leads to transitions between local energy minima (probability maxima).
3. Probabilistic inference using Markov chain Monte Carlo methods by Neal specifically discusses the size of the annealing step.