r语言分别用简单蒙特卡罗估计、对偶蒙特卡罗积分、控制变量蒙特卡罗积分、重要抽样蒙特卡洛积分,分层抽样蒙特卡洛积分估计e^-x/(1+x^2)在[0,1]上的积分并计算方差的代码
时间: 2024-05-20 16:17:51 浏览: 103
# 简单蒙特卡罗估计
set.seed(123)
n <- 10000
x <- runif(n, 0, 1)
y <- exp(-x) / (1 + x^2)
mc_estimate <- mean(y)
mc_var <- var(y) / n
# 对偶蒙特卡罗积分
set.seed(123)
n <- 10000
u <- runif(n, -1, 1)
x <- 1 / (1 + u^2)
y <- exp(-x)
mc_estimate <- mean(y) * 2
mc_var <- var(y) / n
# 控制变量蒙特卡罗积分
set.seed(123)
n <- 10000
x <- runif(n, 0, 1)
y1 <- exp(-x)
y2 <- 1 / (1 + x^2)
cov_xy <- cov(y1, y2)
b <- cov_xy[1, 2] / var(y2)
z <- y1 - b * (y2 - mean(y2))
mc_estimate <- mean(z)
mc_var <- var(z) / n
# 重要抽样蒙特卡洛积分
set.seed(123)
n <- 10000
x <- qcauchy(runif(n))
y <- exp(-x) / (1 + x^2) / dcauchy(x)
mc_estimate <- mean(y)
mc_var <- var(y) / n
# 分层抽样蒙特卡洛积分
set.seed(123)
n <- 10000
m <- 10
h <- 1 / m
x <- runif(n)
idx <- floor(x / h) + 1
y <- rep(0, n)
for (i in 1:m) {
y[idx == i] <- exp(-x[idx == i]) / (1 + x[idx == i]^2) / h
}
mc_estimate <- mean(y)
mc_var <- var(y) / n
阅读全文