Monte Carlo Methods
A general problem in probability and statistical applications is the computation of an expectation of a random variable. We illustrate the use of R to compute some expectations by the Monte Carlo method. These computations are helpful in comparing samplin
- PDF / 474,950 Bytes
- 30 Pages / 439.37 x 666.14 pts Page_size
- 95 Downloads / 238 Views
Monte Carlo Methods
13.1 The Monte Carlo Method of Computing Integrals 13.1.1 Introduction A general problem in probability and statistical applications is the computation of an expectation of a random variable. Suppose the random variable X has a density function f (x) and we are interested in computing the expectation of a function h(x) which is expressible as the integral E(h(X)) = h(x)f (x)dx. In the case where h(x) = x, this expectation corresponds to the mean μ = E(X). If h(x) = (x − μ)2 , then the expectation corresponds to the variance var(X) = E[(X − μ)2 ]. If the function h is the indicator function h(x) = 1 for x ∈ A, and h(x) = 0 elsewhere, then this expectation corresponds to a probability E(h(X)) = P (A). Suppose we are able to simulate a random sample of values x(1) , ..., x(m) from the density f (x). Then the Monte Carlo estimate of the expectation E(h(X)) is given by m h(x(j) ) ¯ = j=1 . h m ¯ will converge to E(h(X)) This is a good estimate of the expectation, since h ¯ is just a sample as the simulation sample size m approaches infinity. But h estimate of the expectation and so there will be some error in this simulation¯ is given by based calculation. The variance of h ¯ = V ar(h(X)) , V ar(h) m
J. Albert and M. Rizzo, R by Example, Use R, DOI 10.1007/978-1-4614-1365-3__13, © Springer Science+Business Media, LLC 2012
307
308
13 Monte Carlo Methods
where V ar(h(X)) is the variance of each simulated draw of h(x(j) ). We estimate V ar(h(X)) by the sample variance of the simulated draws {h(x(j) )}: = V ar(h(X))
m
(j) ¯ 2 j=1 (h(x ) − h)
m−1
.
Then the associated standard error of this Monte Carlo estimate is given by ¯ ≈ V ar(h(X)) . seh¯ = V ar(h) m Example 13.1 (Sleepless in Seattle). Here is a simple example of the Monte Carlo method for computing integrals, inspired by the movie Sleepless in Seattle. Sam and Annie have a rendezvous at the Empire State Building on Valentine’s Day. Sam arrives at a random time between 10 and 11:30 pm and Annie arrives at a random time between 10:30 and 12 am. What is the probability that Annie arrives before Sam? What is the expected difference in the two arrival times?
13.1.2 Estimating a probability Let A and S denote respectively Annie’s and Sam’s arrival times, where we measure the arrival time as the number of hours after noon. We are assuming A and S are independent, where A is uniformly distributed on (10.5, 12) and S is uniformly distributed on (10, 11.5). We wish to compute the probability P (A < S) which is expressed as the integral fA (a)fS (s)dads, P (A < S) = a sam = runif(1000, 10, 11.5) > annie = runif(1000, 10.5, 12)
The probability P (A < S) is estimated by the proportion of simulated pairs (a, s) where a is smaller than s. We count the number of pairs where a < s by the sum function and divide this by the total number of simulations. > prob = sum(annie < sam) / 1000 > prob [1] 0.229
The estimated probability that Annie arrives before Sam is 0.229. Figure 13.1 illustrates the process to compute this probability. The 1000
Data Loading...