掷 n 个骰子后得到特定和的概率。 Ruby

Probability of getting specific sum after rolling n dice. Ruby

求 n 个骰子滚动和的概率的最佳解决方案是什么? 我正在通过找到

来解决它
  1. 平均值。
  2. 标准偏差。
  3. z_score x
  4. 以下的数字
  5. z_score 以上数字 x
  6. 将两者都转换为概率
  7. 相减

这就是我到目前为止所做的。

# sides - number of sides on one die
def get_mean(sides)
  (1..sides).inject(:+) / sides.to_f
end

def get_variance(sides)
  mean_of_squares = ((1..sides).inject {|sum, side| sum + side ** 2}) / sides.to_f
  square_mean = get_mean(sides) ** 2

  mean_of_squares - square_mean
end

def get_sigma(variance)
  variance ** 0.5
end

# x - the number of points in question
def get_z_score(x, mean, sigma)
  (x - mean) / sigma.to_f
end

# Converts z_score to probability
def z_to_probability(z)
  return 0 if z < -6.5
  return 1 if z > 6.5

  fact_k = 1
  sum = 0
  term = 1
  k = 0

  loop_stop = Math.exp(-23)
  while term.abs > loop_stop do
    term = 0.3989422804 * ((-1)**k) * (z**k) / (2*k+1) / (2**k) * (z**(k+1)) / fact_k
    sum += term
    k += 1
    fact_k *= k
  end

  sum += 0.5
  1 - sum
end

# Calculate probability of getting 'х' total points by rolling 'n' dice with 'sides' number of sides.
def probability_of_sum(x, n, sides=6)

  mean = n * get_mean(sides)
  variance = get_variance(sides)
  sigma = get_sigma(n * variance)

  # Rolling below the sum
  z1 = get_z_score(x, mean, sigma)
  prob_1 = z_to_probability(z1)

  # Rolling above the sum
  z2 = get_z_score(x+1, mean, sigma)
  prob_2 = z_to_probability(z2)

  prob_1 - prob_2
end

# Run probability for 100 dice
puts probability_of_sum(400, 100)

让我担心的是,当我选择x = 200时,概率为0。 这是正确的解决方案吗?

我也是用Monte Carlo的方法解决的,结果比较接近

# x - sum of points to find probability for
# n - number of dice
# trials - number of trials
def monte_carlo(x, n, trials=10000)
  pos = 0

  trials.times do
    sum = n.times.inject(0) { |sum| sum + rand(1..6) }
    pos += 1 if sum == x
  end

  pos / trials.to_f
end

puts monte_carlo(300, 100, 30000)

将两个独立概率分布的结果相加与convolving两个分布的结果相同。如果分布是离散的,那么它是一个离散卷积。

因此,如果单个骰子表示为:

probs_1d6 = Array.new(6) { Rational(1,6) }

那么2d6可以这样计算:

probs_2d6 = []
probs_1d6.each_with_index do |prob_a,i_a|  
  probs_1d6.each_with_index do |prob_b,i_b| 
    probs_2d6[i_a + i_b] = ( probs_2d6[i_a + i_b] || 0 ) + prob_a * prob_b
  end
end

probs_2d6
# =>  [(1/36), (1/18), (1/12), (1/9), (5/36), (1/6), 
#      (5/36), (1/9), (1/12), (1/18), (1/36)]

虽然这是骰子边的 n 平方,并且完全合乎逻辑的组合可以减少它,但对于更复杂的设置,这样做通常不太灵活。这种方法的好处是您可以继续添加更多骰子并进行其他更奇特的组合。例如,要获得 4d6,您可以对 2d6 的两个结果进行卷积。使用有理数可以避免浮点精度问题。

我跳过了一个细节,您确实需要存储初始偏移量(对于普通的六面骰子为 +1)并将其加在一起,以便了解概率匹配。

我在 gem games_dice 中用浮点而不是 Rational 制作了这个逻辑的更复杂的版本,它可以处理其他一些骰子组合。

这是使用上述方法以一种天真的方式对您的方法进行的基本重写(一次简单地组合一个骰子的效果):

def probability_of_sum(x, n, sides=6)
  return 0 if x < n
  single_die_probs = Array.new(sides) { Rational(1,sides) }

  combined_probs = [1] # Represents chance of getting 0 when rolling 0 dice :-)

  # This is not the most efficient way to do this, but easier to understand
  n.times do
    start_probs = combined_probs
    combined_probs = []
    start_probs.each_with_index do |prob_a,i_a|  
        single_die_probs.each_with_index do |prob_b,i_b| 
          combined_probs[i_a + i_b] = ( combined_probs[i_a + i_b] || 0 ) + prob_a * prob_b
        end
    end
  end

  combined_probs[ x - n ] || 0
end

puts probability_of_sum(400, 100).to_f
# => 0.0003172139126369326

注意该方法实际上计算的是100到600的全概率分布,所以你只需要调用一次,存储一次数组(加上偏移量+100),就可以做概率等其他有用的事情了变得大于一定数量。由于在 Ruby.

中使用了 Rational 个数字,所有这些都具有完美的精度

因为在你的情况下,你只有一种骰子,我们可以避免使用 Rational 直到最后,只使用整数(基本上是组合值的计数),然后除以总数组合(边数为卷数的幂)。这要快得多,并且 returns 不到一秒便可获得 100 个骰子的值:

def probability_of_sum(x, n, sides=6)
  return 0 if x < n
  combined_probs = [1] # Represents chance of getting 0 when rolling 0 dice :-)

  n.times do
    start_probs = combined_probs
    combined_probs = []
    start_probs.each_with_index do |prob_a,i_a|  
        sides.times do |i_b| 
          combined_probs[i_a + i_b] = ( combined_probs[i_a + i_b] || 0 ) + prob_a
        end
    end
  end

  Rational( combined_probs[ x - n ] || 0, sides ** n )
end

存在一个涉及二项式系数交替求和的精确解。我已经在几个地方(在 Quora and MSE 上)写了它,你可以在别处找到它,尽管有一些有缺陷的版本。请注意,如果你实现它,你可能需要取消比最终结果大得多的二项式系数,如果你使用浮点运算,你可能会失去太多的精度。

Neil Slater 关于使用动态规划计算卷积的建议是一个很好的建议。它比二项式系数的求和慢,但相当稳健。您可以通过几种方式加快速度,例如通过平方求幂,以及使用快速傅里叶变换,但很多人会发现这些方法有些矫枉过正。

要修正你的方法,你应该使用(简单的)连续性修正法线近似,并限制在你有足够骰子的情况下,你的评估距离你期望的法线近似值的最大值和最小值足够远近似值将是好的,无论是在绝对意义上还是相对意义上。连续性修正就是你用n-1/2到n+1/2的区间代替n的计数。

总共掷出200的方法数的精确计数是7745954278770349845682110174816333221135826585518841002760,所以概率是除以6^100,大约是1.18563 x 10^-20。

具有简单连续性校正的正态近似为 Phi((200.5-350)/sqrt(3500/12))-Phi((199.5-350)/sqrt(3500/12)) = 4.2 x 10^ -19。这在绝对值上是准确的,因为它非常接近于 0,但它偏离了 35 倍,因此相对而言并不是很好。正态近似给出了更接近中心的更好的相对近似。

这是我的最终版本。

  1. get_z_score 中的总和偏移量分别更改为 x-0.5x+0.5 以获得更精确的结果。
  2. 添加 return 0 if x < n || x > n * sides 以涵盖总和小于骰子数且大于骰子数乘以面数的情况。
  3. 添加了带有结果的基准测试

主要功能

# sides - number of sides on one die
def get_mean(sides)
  (1..sides).inject(:+) / sides.to_f
end

def get_variance(sides)
  mean_of_squares = ((1..sides).inject {|sum, side| sum + side ** 2}) / sides.to_f
  square_mean = get_mean(sides) ** 2

  mean_of_squares - square_mean
end

def get_sigma(variance)
  variance ** 0.5
end

# x - the number of points in question
def get_z_score(x, mean, sigma)
  (x - mean) / sigma.to_f
end

# Converts z_score to probability
def z_to_probability(z)
  return 0 if z < -6.5
  return 1 if z > 6.5

  fact_k = 1
  sum = 0
  term = 1
  k = 0

  loop_stop = Math.exp(-23)
  while term.abs > loop_stop do
    term = 0.3989422804 * ((-1)**k) * (z**k) / (2*k+1) / (2**k) * (z**(k+1)) / fact_k
    sum += term
    k += 1
    fact_k *= k
  end

  sum += 0.5
  1 - sum
end

# Calculate probability of getting 'х' total points by rolling 'n' dice with 'sides' number of sides.
def probability_of_sum(x, n, sides=6)
  return 0 if x < n || x > n * sides

  mean = n * get_mean(sides)
  variance = get_variance(sides)
  sigma = get_sigma(n * variance)

  # Rolling below the sum
  z1 = get_z_score(x-0.5, mean, sigma)
  prob_1 = z_to_probability(z1)

  # Rolling above the sum
  z2 = get_z_score(x+0.5, mean, sigma)
  prob_2 = z_to_probability(z2)

  prob_1 - prob_2
end

基准

require 'benchmark'

Benchmark.bm do |x|
  x.report { @prob = probability_of_sum(350, 100).to_f }
  puts "\tWith x = 350 and n = 100:"
  puts "\tProbability: #{@prob}"
end
puts

Benchmark.bm do |x|
  x.report { @prob = probability_of_sum(400, 100).to_f }
  puts "\tWith x = 400 and n = 100:"
  puts "\tProbability: #{@prob}"
end
puts

Benchmark.bm do |x|
  x.report { @prob = probability_of_sum(1000, 300).to_f }
  puts "\tWith x = 1000 and n = 300:"
  puts "\tProbability: #{@prob}"
end

结果

       user     system      total        real
   0.000000   0.000000   0.000000 (  0.000049)
    With x = 350 and n = 100:
    Probability: 0.023356331366255034

       user     system      total        real
   0.000000   0.000000   0.000000 (  0.000049)
    With x = 400 and n = 100:
    Probability: 0.00032186531055478085

       user     system      total        real
   0.000000   0.000000   0.000000 (  0.000032)
    With x = 1000 and n = 300:
    Probability: 0.003232390001131513