MCMC Practice

STA602 at Duke University

Metropolis algorithm

Exercise

Suppose the target distribution we wish to sample from is given by probability mass function

\[ \pi(\theta) = \theta / w \text{ for } \theta \in \{1, 2, \ldots 6\} \]

in words, we wish to roll a die with probability \(1/w\) of landing on face 1, \(2/w\) of landing on face 2, etc.

  • Write a Metropolis algorithm to approximate the target distribution using a proposal \(J(\theta = j | \theta^{(s)} = i) = 1/6\) for all \(j\), i.e. propose a new state \(j\) uniformly. Run your Markov chain for \(S=10000\) states.

  • The Metropolis algorithm requires a symmetric proposal \(J\). Explain why this proposal is symmetric.

  • Plot a histogram of the Markov chain samples. Does the plot match your intuition?

  • Compare the estimated probabilities of each outcome to the truth (compute \(w\)).

Metropolis-Hastings

Metropolis-Hastings algorithm

Let \(\pi(\theta)\) be the target distribution. The Metropolis-Hastings algorithm proceeds:

  1. sample \(\theta^{*} \sim J(\theta | \theta^{(s)})\);

  2. compute the acceptance ratio

\[ r = \frac{\pi(\theta^*)}{\pi(\theta^{(s)})} \times \frac{J(\theta^{(s)}| \theta^*)}{ J(\theta^{*}| \theta^{(s)}) } \]

  1. set \(\theta^{(s+1)}\) to \(\theta^*\) with probability \(\min(1, r)\), otherwise set \(\theta^{(s+1)}\) to \(\theta^{(s)}\).

Important: We correct for asymmetry; the proposal distribution \(J\) need not be symmetric!

Exercise

Metropolis-Hastings lets us work with non-symmetric proposals. Re-write the algorithm of the previous exercise using the non-symmetric proposal \(J(\theta = j | \theta^{(s)} = i)\) such that

\[ \theta = \begin{cases} 1 & \text{ with prob } & 0.05\\ 2 & \text{ with prob } & 0.15\\ 3 & \text{ with prob } & 0.2\\ 4 & \text{ with prob } & 0.15\\ 5 & \text{ with prob } & 0.15\\ 6 & \text{ with prob } & 0.3\\ \end{cases} \]

  • compare your results to that those of the previous exercise. In particular, compare the ESS of \(\theta\) in each chain. Which do you prefer? How might you explain this difference in ESS?