r - 蒙特卡罗模拟 - 怎么了？

Question

我想使用函数的蒙特卡罗模拟计算曲线下的面积

. 我想在区间 [-2, 2] 上计算它

到目前为止我的工作

# Define function f
f <- function(x) x^2 + 1

# I want to close my function in rectangle (-2, 2) - x axis and (1, 5) y -axis
n <- 10^6
# Randomize from x axis of rectangle
x_n <- runif(n, min = -2, max = 2)
# Randomize from y axis of rectangle
y_n <- runif(n, min = 1, max = 5)
# Calculate function values of randomized points
values <- f(x_n)

# Formula for are under the curve for monte carlo simulation is 
# Area of rectangle * (Points below curve) / (Total number of points)

所以我的结果是：

> sum(y_n < values) / n * (4 * 4)
[1] 5.329888

这是不好的结果（正确的结果是 9.33333）。我做错了什么？当然，经过百万次采样后，算法应该更接近 9.3333

标签： rsimulationmontecarlointegral

Answer 1

这是一个显示您正在使用的内容的情节。我希望它能帮助您更好地理解我在评论中写的内容：

您似乎忽略了 y=1 下方的矩形。它的面积 (=4) 是缺失的数量。所以代码对于计算非偏移表达式 x^2 是正确的。更改为 y_n <- runif(n, min = 0, max = 5) 并重新运行计算

评论是答案的一半，即您没有模拟 y_n 的 0 和 1 之间的点。那些需要在蒙特卡洛模型中整合一个区域。另一个修改是将 [-2 < x <2]x[0<y<5] = 4*5 的正确总面积添加到正在考虑的总“面积”的计算中

f <- function(x) x^2 + 1

# I want to close my function in rectangle (-2, 2) - x axis and (1, 5) y -axis
n <- 10^6
# Randomize from x axis of rectangle
x_n <- runif(n, min = -2, max = 2)
# Randomize from y axis of rectangle
y_n <- runif(n, min = 0, max = 5)
# Calculate function values of randomized points
values <- f(x_n)

# Formula for are under the curve for monte carlo simulation is 
# Area of rectangle * (Points below curve) / (Total number of points)
 sum(y_n < values) / n * (5 * 4)
#[1] 9.3429  inaccurate to 1 or 2 parts in 933

显示第二种情况的 100 点图：

您可能考虑的另一个 mod 是使用 set.seed 使您的计算可重现。