Acceptance-rejection methods

This post is based on chapter 1.4 of Advanced Markov Chain Monte Carlo.

Another method of generating random variates from distributions is to use acceptance-rejection methods. Basically to generate a random number from , we generate a RN from an envelope distribution , where .  The acceptance-rejection algorithm is as follows:

Repeat until we generate a value from step 2:

1. Generate from and from

2. If , return (as a random deviate from ).

Example: the standard normal distribution

This example illustrates how we generate RNs using the logistic distribution as an envelope distribution. First, note that

On setting , we get . This method is fairly efficient and has an acceptance rate of

since both and are normalised densities.

R code

This example is straightforward to code:

myrnorm = function(M){
  while(1){
    u = runif(1); x = rlogis(1, scale = 0.648)
    if(u < dnorm(x)/M/dlogis(x, scale = 0.648))
      return(x)
  }
}

To check the results, we could call myrnorm a few thousand times:

hist(replicate(10000, myrnorm(1.1)), freq=FALSE)
lines(seq(-3, 3, 0.01), dnorm(seq(-3, 3, 0.01)), col=2)

Example: the standard normal distribution with a squeeze

Suppose the density is expensive to evaluate. In this scenario we can employ an easy to compute function , where . is called a squeeze function. In this example, we’ll use a simple rectangular function, where for . This is shown in the following figure:

The modified algorithm is as follows:

Repeat until we generate a value from step 2:

1. Generate from and from

2. If or , return (as a random deviate from ).

Hence, when we don’t have to compute . Obviously, in this example isn’t that difficult to compute.