Monte Carlo 使用给定提案函数的重要性采样进行集成
Monte Carlo integration using importance sampling given a proposal function
给定拉普拉斯分布方案:
g(x) = 1/2*e^(-|x|)
和样本大小 n = 1000
,我想进行 Monte Carlo (MC) 积分以估计 θ:
通过重要性抽样。最后我想在 R 中计算这个 MC 估计值的均值和标准差。
编辑(在下面的答案之后迟到了)
到目前为止,这是我的 R 代码:
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)
现在我们可以编写一个简单的 R 函数来从拉普拉斯分布中采样:
## `n` is sample size
rlaplace <- function (n) {
u <- runif(n, 0, 1)
ifelse(u < 0.5, log(2 * u), -log(2* (1 - u)))
}
也写一个拉普拉斯分布的密度函数:
g <- function (x) ifelse(x < 0, 0.5 * exp(x), 0.5 * exp(-x))
现在,您的被积函数是:
f <- function (x) {
ifelse(x > 0, exp(-sqrt(x) - 0.5 * x) * sin(x) ^ 2, 0)
}
现在我们使用 1000 个样本估计积分(set.seed
用于再现性):
set.seed(0)
x <- rlaplace(1000)
mean(f(x) / g(x))
# [1] 0.2648853
也与使用积分的数值积分进行比较:
integrate(f, lower = 0, upper = Inf)
# 0.2617744 with absolute error < 1.6e-05
给定拉普拉斯分布方案:
g(x) = 1/2*e^(-|x|)
和样本大小 n = 1000
,我想进行 Monte Carlo (MC) 积分以估计 θ:
通过重要性抽样。最后我想在 R 中计算这个 MC 估计值的均值和标准差。
编辑(在下面的答案之后迟到了)
到目前为止,这是我的 R 代码:
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)
现在我们可以编写一个简单的 R 函数来从拉普拉斯分布中采样:
## `n` is sample size
rlaplace <- function (n) {
u <- runif(n, 0, 1)
ifelse(u < 0.5, log(2 * u), -log(2* (1 - u)))
}
也写一个拉普拉斯分布的密度函数:
g <- function (x) ifelse(x < 0, 0.5 * exp(x), 0.5 * exp(-x))
现在,您的被积函数是:
f <- function (x) {
ifelse(x > 0, exp(-sqrt(x) - 0.5 * x) * sin(x) ^ 2, 0)
}
现在我们使用 1000 个样本估计积分(set.seed
用于再现性):
set.seed(0)
x <- rlaplace(1000)
mean(f(x) / g(x))
# [1] 0.2648853
也与使用积分的数值积分进行比较:
integrate(f, lower = 0, upper = Inf)
# 0.2617744 with absolute error < 1.6e-05