Fitting a Mixture Sampling Model, Part II

In the previous post, we introduced the mixture sampling model and showed that the posterior had an interesting bimodal shape. Here we illustrate the use of two Metropolis random walk algorithms to sample from this posterior.

Recall that the definition of the log posterior is in poisson.mix.R and the data is contained in the vector y.

We begin with a random walk algorithm using an identity var-cov matrix and a scale constant of 0.01. We begin at value (3, 3) and run for 10,000 iterations. We store the fit in the variable fit1.

proposal=list(var=diag(c(1,1)),scale=.01)
start=array(c(3,3),c(1,2))
fit1=rwmetrop(poisson.mix,proposal,start,10000,list(y=y,p=.4))

In our second algorithm, we also use an identity var-cov matrix, but use the large scale constant of 0.2.

proposal=list(var=diag(c(1,1)),scale=.2)
start=array(c(3,3),c(1,2))
fit2=rwmetrop(poisson.mix,proposal,start,10000,list(y=y,p=.4))

The acceptance rates of the two algorithms are:

Acceptance rate
fit1 94%
fit2 23%

One basic method of MCMC output analysis uses trace plots which are time series plots of the simulated draws. We show trace plots for each random walk below.

Random walk (scale = 0.01)

Note the snake-like appearance of this graph. This is typical of MCMC random walk runs with high acceptance rates.

Random walk (scale = 0.2)

This plot displays better mixing and it visits both areas of concentration of the posterior.

To check the accuracy of each chain in simulating the marginal posterior of theta1 = log(lambda1), we first use the function simcontour that is based on simulations from the grid that covers the actual posterior. A density estimate on the draws of theta1 based on this "exact" algorithm is presented as a red line and the density estimate of the draws of theta1 from the MCMC algorithm are presented in blue.

We first show the random walk with scale 0.01 and then the random walk with scale 0.2.



In this example, it seems that the second MCMC run produced a better approximation to the marginal posterior density of theta1. The first run with scale 0.01 actually missed the important area of the posterior.