R Programming - 掷骰子算法答案

【问题标题】：R Programming - throwing dice algorithmR Programming - 掷骰子算法
【发布时间】：2017-08-24 16:14:52
【问题描述】：

我想在 R 中运行一个程序，要求用户选择一些骰子，在骰子上运行模拟并确定将最小数字滚动到最大数字的概率。

例如，如果用户选择 5 个骰子，则最小掷骰数为 5x1=5，最大掷骰数为 5x6=30。我已经有一组骰子和一组总数的代码 - 只需要知道如何增加它。 “d”是骰子数，“k”是掷骰总数，“nreps”是模拟运行（例如 1,000,000）。我想将每个概率存储在一个向量中，然后给出概率与滚动总数（从最小值到最大值）的图（泊松分布）。

probtotk <- function(d, k, nreps){
  count <- 0
  #do the experiment nreps times
  for (rep in 1:nreps){
    total <- sum(sample(1:6, d, replace = TRUE))
  if (total == k) count <- count +1
}
 return(count/nreps) 
}

【问题讨论】：

我对你的目标感到困惑。您似乎想关注一个特定的总数（k），但也想查看每个可能的总数（“给出概率与滚动总数的图（从最小值到最大值）”我>）。似乎您应该只选择这两个选项之一（从最小值到最大值似乎更有趣，也同样简单）。
我的意思是，如果你掷 5 个骰子，那么进行 1M 模拟以找到总数为 5 的概率，然后进行 1M 模拟以找到总数为 6 的概率是没有意义的, ... 30. 相反，只需进行 1M 次模拟并查看总数的分布。

标签： r

【解决方案1】：

我们可以使用 R 的矢量化来非常快速地做到这一点。正如我的 cmets 所建议的，我不会使用k。

对于d 骰子和nreps 模拟，我们将有d * nreps 的总骰数。我们使用sample(6, size = d * nreps, replace = T) 一次性模拟所有这些。我们将结果放在一个矩阵中，该矩阵包含nreps 列和d 行，因此每一列代表一卷d 骰子。列总和给出每卷的总数。 table函数统计每个总数的出现次数，prop.table函数将其转化为比例。

dice_tot_prob = function(d, nreps) {
    rolls = matrix(sample(6, size = d * nreps, replace = T), ncol = nreps)
    totals = colSums(rolls)
    return(prop.table(table(totals)))
}

dice_tot_prob(5, 1e5)
totals
      5       6       7       8       9      10      11      12      13      14      15      16      17      18      19      20      21      22      23 
0.00015 0.00066 0.00200 0.00446 0.00904 0.01615 0.02655 0.03958 0.05456 0.07013 0.08379 0.09511 0.10065 0.10068 0.09214 0.08391 0.06936 0.05384 0.03891 
     24      25      26      27      28      29      30 
0.02576 0.01664 0.00880 0.00474 0.00180 0.00044 0.00015

prop.table 结果很好，因为它有一个默认的绘图方法：

plot(dice_tot_prob(5, 1e5))

【讨论】：

相同的解决方案但更好的解释:-)
嘿 Gregor...我才刚接触 R 语言几天。这太酷了！我可以在这里说“怪胎”吗？谢谢！
输入了您的代码，但找不到函数“prob.table”。我需要下载一个包吗？
复制/粘贴比打字更安全。它是prop.table 和p，如proportion。
键入它以了解所有这些东西是如何编写的。进行了修复，并且可以正常工作。感谢您的所有帮助。

【解决方案2】：

我想你需要的是这样的：

library(magrittr)
sample(1:6, nreps * d, replace = TRUE) %>%
  matrix(nrow = d) %>%
  colSums() %>%
  table() %>%
  divide_by(nreps)

【讨论】：