Search code examples
rsamplingmontecarlointegralnumerical-integration

Monte Carlo integration using importance sampling given a proposal function


Given a Laplace Distribution proposal:

g(x) = 1/2*e^(-|x|)

and sample size n = 1000, I want to Conduct the Monte Carlo (MC) integration for estimating θ:

enter image description here

via importance sampling. Eventually I want to calculate the mean and standard deviation of this MC estimate in R once I get there.


Edit (arrived late after the answer below)

This is what I have for my R code so far:

library(VGAM)
n = 1000
x = rexp(n,0.5)
hx = mean(2*exp(-sqrt(x))*(sin(x))^2)
gx = rlaplace(n, location = 0, scale = 1)

Solution

  • enter image description here

    Now we can write a simple R function to sample from Laplace distribution:

    ## `n` is sample size
    rlaplace <- function (n) {
      u <- runif(n, 0, 1)
      ifelse(u < 0.5, log(2 * u), -log(2* (1 - u)))
      }
    

    Also write a function for density of Laplace distribution:

    g <- function (x) ifelse(x < 0, 0.5 * exp(x), 0.5 * exp(-x))
    

    Now, your integrand is:

    f <- function (x) {
      ifelse(x > 0, exp(-sqrt(x) - 0.5 * x) * sin(x) ^ 2, 0)
      }
    

    Now we estimate the integral using 1000 samples (set.seed for reproducibility):

    set.seed(0)
    x <- rlaplace(1000)
    mean(f(x) / g(x))
    # [1] 0.2648853
    

    Also compare with numerical integration using quadrature:

    integrate(f, lower = 0, upper = Inf)
    # 0.2617744 with absolute error < 1.6e-05