R语言 蒙特卡洛模拟计算π

2020-12-29 15:29:09 浏览数 (1)

蒙特卡洛法的基本思想是：为了求解问题，首先建立一个概率模型或随机过程，使它的参数或数字特征等于问题的解：然后通过对模型或过程的观察或抽样试验来计算这些参数或数字特征，最后给出所求解的近似值。解的精确度用估计值的标准误差来表示。蒙特卡洛法的主要理论基础是概率统计理论，主要手段是随机抽样、统计试验。用蒙特卡洛法求解实际问题的基本步骤为：

  1. 根据实际问题的特点．构造简单而又便于实现的概率统计模型．使所求的解恰好是所求问题的概率分布或数学期望；
  2. 给出模型中各种不同分布随机变量的抽样方法；
  3. 统计处理模拟结果，给出问题解的统计估计值和精度估计值。

考虑平面上的一个边长为1的正方形及其内部的一个形状不规则的“图形”，如何求出这个“图形”的面积呢？Monte Carlo方法是这样一种“随机化”的方法：向该正方形“随机地”投掷N个点，有M个点落于“图形”内，则该“图形”的面积近似为M/N。

代码语言：javascript复制

center <- c(1,1)
distanceToCenter <- function(a){sqrt(sum((center-a)^2))}

set.seed(1234)
n <- 10000 #设定生成n个点
mx <- matrix(runif(n*2,min = 0,max = 2),nrow = n,ncol = 2,byrow = T)
b <- apply(mx,1,distanceToCenter)
# 求向量distance中元素值小于半径1的个数
# 因逻辑表达式结果要么0，要么1，故可用distance<1
num <- mean(b<1)
pi <= 4*num
pi

par(bg="beige")
plot(mx,col="azure3",asp = 1)
abline(h=0,col="red",lty="dotdash",lwd=2)
abline(h=2,col="red",lty="dotdash",lwd=2)
abline(v=0,col="red",lty="dotdash",lwd=2)
abline(v=2,col="red",lty="dotdash",lwd=2)


points(mx[b<1,],col="green")

library(plotrix)
draw.circle(1,1,1,border="coral2",lty="dashed",lwd=2) #绘制一个圆
points(x=1,y=1,col="red",pch=20,cex=1.5,lwd=1.5)
统计

0 人点赞