Robert_V._Hogg,_Joseph_W._McKean,_Allen_T._Craig

(Jacob Rumans) #1
4.8. The Method of Monte Carlo 299

Hence, by differentiating both sides, we find that the pdf ofXisf(x).

There are two facts worth noting. First, the probability of an acceptance in the
algorithm is 1/M. This can be seen in the derivation in the proof of the theorem.
Just consider the denominators in the derivation which show that


P

[
U≤
f(Y)
Mg(Y)

]
=
1
M

. (4.8.8)


Hence, for efficiency of the algorithm we wantMas small as possible. Secondly,
normalizing constants of the two pdfsf(x)andg(x) can be ignored. For example,
iff(x)=kh(x)andg(x)=ct(x)forconstantscandk, then we can use the rule


h(x)≤M 2 t(x), −∞<x<∞, (4.8.9)

and change the ratio in step 2 of the algorithm toU≤h(Y)/[M 2 t(Y)]. It follows
directly that expression (4.8.5) holds if and only if expression (4.8.9) holds where
M 2 =cM/k. This often simplifies the use of the accept–reject algorithm.
We next present two examples of the accept–reject algorithm. The first exam-
ple offers a normal generator where the instrumental random variable,Y,hasa
Cauchy distribution. The second example shows how all gamma distributions can
be generated.


Example 4.8.7.Suppose thatXis a normally distributed random variable with
pdfφ(x)=(2π)−^1 /^2 exp{−x^2 / 2 }andYhas a Cauchy distribution with pdfg(x)=
π−^1 (1 +x^2 )−^1. As Exercise 4.8.9 shows, the Cauchy distribution is easy to simulate
because its inverse cdf is a known function. Ignoring normalizing constants, the
ratio to bound is


f(x)
g(x)

∝(1 +x^2 )exp{−x^2 / 2 }, −∞<x<∞. (4.8.10)

As Exercise 4.8.17 shows, the derivative of this ratio is−xexp{−x^2 / 2 }(x^2 −1),
which has critical values at±1. These values provide maxima to (4.8.10). Hence,


(1 +x^2 )exp{−x^2 / 2 }≤2exp{− 1 / 2 }=1. 213 ,

soM 2 =1.213. Hence, from the above discussion,M=(π/


2 π)1.213 = 1. 520.
Hence, the acceptance rate of the algorithm is 1/M=0.6577.

Example 4.8.8.Suppose we want to generate observations from a Γ(α, β). First,
ifYhas a Γ(α,1)-distribution thenβY has a Γ(α, β)-distribution. Hence, we need
only consider Γ(α,1) distributions. So letXhave a Γ(α,1)-distribution. Ifαis a
positive integer then by Theorem 3.3.1 we can writeXas

X=T 1 +T 2 +···+Tα,

whereT 1 ,T 2 ,···,Tαare independent and identically distributed with the common
Γ(1,1)-distribution. In the discussion around expression (4.8.2), we have shown how
to generateTi.

Free download pdf