Random walk methods: the Metropolis algorithm

Suppose that we want to generate random variables according to an arbitrary probability density $P(x)$. The Metropolis algorithm produces a ``random walk'' of points $\{x_i\}$ whose asymptotic probability approaches $P(x)$ after a large number of steps. The random walk is defined by a ``transition probability'' $w(x_i \rightarrow x_j)$ for one value $x_i$ to another $x_j$ in order that the distribution of points $x_0$, $x_1$, $x_2$, ... converges to $P(x)$. In can be shown that it is sufficient (but not necessary) to satisfy the ``detailed balance'' condition

$$p(x_i)w(x_i \rightarrow x_j) = p(x_j)w(x_j \rightarrow x_i).$$ (265)

This relation dos not specify $w(x_i \rightarrow x_j$ uniquely. A simple choice is

$$w(x_i \rightarrow x_j)=\min{\left[ 1,\frac{P(x_j)}{P(x_i)} \right] }.$$ (266)

This choice can be described by the following steps. Suppose that the ``random walker'' is a position $x_n$. To generate $x_{n+1}$ we

1. choose a trial position $x_t=x_n+\delta _n$ , where the $\delta _n$ is a random number in the interval $[-\delta ,\delta]$.

2. Calculate $w=P(x_t)/P(x_n)$.

3. If $w \geq 1$ we accept the change and let $x_{n+1}=x_t$.

4. If $w \leq 1$, generate a random number $r$.

5. If $r \leq w$, accept the change and let $x_{n+1}=x_t$.

6. If the trial change is not accepted, the let $x_{n+1}=x_n$.

It is necessary to sample a number of points of the random walk before the asymptotic probability $P(x)$ is attained. How do we choose the ``step size'' $\delta$? If $\delta$ is too large, only a small fraction of changes will be accepted and the sampling will be inefficient. If $\delta$ is too small, a large number will be accepted, but it would take too long to sample $P(x)$ over the whole interval of interest. Ideally, we want at least 1/3-1/2 of the trial steps to be accepted. We also want to choose $x_0$ such that the distribution $\{x_i\}$ converges to $P(x)$ as quickly as possible. An obvious choice is to begin the random walk at the point where $P(x)$ is maximum.