R中的编程概率分布

机器算法验证 r 可能性 自习 分布
2022-03-24 04:20:53

我试图在 R 中模拟这个问题的答案。

明尼苏达州西北部春季容易发生洪水。假设 20% 的农民投保了洪水损失保险。随机选择四名农民。至少有两个农民拥有洪水保险的概率是多少?

感谢任何人的帮助。

2个回答

这个答案将引导您完成大多数随机模拟常见的步骤,展示如何在任何软件平台上以小、简单、易于测试的块创建代码。该过程用R代码说明。您可以一边运行代码片段,一边查看它们产生的结果。

FWIW,这是一种紧凑、快速且肮脏的解决方案,可将计算压缩成一行。它打印出基于n迭代模拟的估计值。(计算 CPU 时间每秒大约 400 万次迭代。)

n <- 1e6; mean(rowSums(matrix(runif(4*n) <= 0.20, nrow=n)) >= 2)

第 1 步:模拟一位农民。

创建一个带有门票的盒子,每个农民一个。那些票上写着“已投保”和“未投保”。盒子里“被保险人”的比例必须是20%。(有关随机变量的盒装票模型的详细说明,请参阅随机变量的含义

在几乎所有支持随机数生成的软件系统中,都有一个功能可以从一个非常大的盒子中抽出一张票(并在看到票后更换它)。这些票的浮点值从到略小于对于任何区间,值在之间的票的比例是之间的所有数字的 20% 上写“保险”来利用此功能例如,您可以在之间的所有数字上写上“保险” 。010ab<1abba01020/100

这是均匀分布(在之间)。这个函数中被调用从这个盒子中抽一次票可以被编程为01Rrunif

label <- ifelse(runif(1) <= 20/100, "insured", "not insured")

但是,取消标签更快更方便,并将其替换为数字 (为所有其他票保留这给出了更简单的代码10R

label.indicator <- runif(1) <= 20/100

因为在R(在许多系统中)一个真值等于和一个假值10

第 2 步:模拟四个农民。

只需独立从盒子里抽票。这是通过循环完成的。R您从runif. 重要的是要知道这些值是(伪)独立的:它们似乎并不相互依赖。因此,

label.indicators <- runif(4) <= 20/100

将模拟从这个盒子中抽取四个农民的票。FALSE它从集合 { , TRUE} 或等效中生成一个由四个数字组成的数组{0,1}

第 3 步:计算统计数据(农民人数)。

任何组中的农民人数都是通过将他们的指标相加得出的。

farmer.count <- sum(label.indicators)

这是因为每个有保险的农民为总和贡献了,每个没有保险的农民贡献了该金额仅计入参保农民。使用指标而不是标签的效率在这里是显而易见的。10

第四步:重复多次。

这是大多数平台中的循环。在某些情况下,包括R,抽出很多票并将它们分成四组会更快。(这是因为生成许多随机值的速度通常几乎与生成一个随机值一样快:所涉及的开销更少。)每个组代表一次模拟迭代。以下代码将每次迭代的四张票放入数组的行中,然后对每一行求和(如步骤 3 中所示)。

n.sim <- 1e6   # Number of iterations
n.farmers <- 4 # Number of farmers per iteration
simulation <- rowSums(matrix(runif(n.sim * n.farmers) <= 20/100, nrow=n.sim))

第 5 步:后处理。

观察到两个或更多被保险农民的机会可以估计为在样本中发现两个或更多被保险农民的迭代比例。和以前一样,这可以通过测试和求和(或测试和平均)有效地找到:

estimate <- mean(simulation >= 2)

第 6 步:评估。

您刚刚进行了计算机实验。与任何其他实验一样,结果是可变的,因此您应该为结果提供标准误差。在这个二项式实验中,估计值(有次迭代)的标准误差是p^n

se(p^)=p^(1p^)/n.

计算并打印结果:

se.estimate <- sqrt(estimate * (1-estimate) / n.sim)
print(c(Estimate=estimate, SE=se.estimate), digits=2)

当随机种子设置为(请参阅下面的完整代码)并运行一百万次迭代时,输出为17

Estimate       SE 
 0.18127  0.00039 

这是完整的解决方案。 它的结构允许参数的轻松变化,因此,通过重复它(计算时间不到一秒),您可以研究结果如何依赖于参数,从而更直观地了解正在发生的事情。

#
# Specify the problem.
#
p <- 20/100    # Chance of insurance
k <- 2         # Minimum number of insured farmers
n.farmers <- 4 # Number of farmers per iteration
n.sim <- 1e6   # Number of iterations
#
# Simulate.
#
set.seed(17)   # Optional: provides a reproducible result
simulation <- rowSums(matrix(runif(n.sim * n.farmers) <= p, nrow=n.sim))
#
# Post-process and report.
#
estimate <- mean(simulation >= k)
se.estimate <- sqrt(estimate * (1-estimate) / n.sim)
print(c(Estimate=estimate, SE=se.estimate), digits=2)

这是您问题的可能解决方案

prop=0
T=1e6
for (t in 1:T){
  insrd=sample(0:1,4,rep=TRUE,prob=c(8,2))
  prop=prop+(sum(insrd)>1)}
print(prop/T)

答案为 0.1809。但是,当计算四名投保农民的数量时,您可以通过意识到四名农民的抽签是二项式 B(4,0.2) 随机变量来计算该概率的确切值 [0.1808]。