Multi-target guessing probability.

multitarget.rb

# This script answers the following question:
# Alice chooses N random numbers between 1 and K.
# Bob chooses G random numbers between 1 and K.
# What is the probability that at least one number is chosen by both of them?

# Computes (K-N choose G) / (K choose G) in O(N)-ish time.

k = 1_000_000_000
n = 10_000
g = 100_000

div = 1
mul = 1

n.times do |i|
  div *= (k - i)
  mul *= (k - g - i)
end

# Why does this work?
# If you work it out...
#   (K-N choose G) / (K choose G)
# ...is equivalent to...
#   (K-N)!(K-M)! / ( K! (K-N-M)! )
# ...which can be understood like this (nums in box are the exponent):
#
#  [      0       |      0       |       0      |     0    ]
#  K             K-N            K-M           K-N-M        0
#
# Apply the K! in the denominator, and it becomes:
#
#  [      -1      |      -1      |      -1      |    -1    ]
#  K             K-N            K-M           K-N-M        0
# 
# Apply the (K-N-M)! in the denominator, and it becomes:
#
#  [      -1      |      -1      |      -1      |    -2    ]
#  K             K-N            K-M           K-N-M        0
#
# Apply the (K-N)! in the numerator and it becomes:
#
#  [      -1      |       0      |       0      |    -1    ]
#  K             K-N            K-M           K-N-M        0
#
# Apply the (K-M)! in the numerator and it becomes:
#
#  [      -1      |       0      |       1      |     0    ]
#  K             K-N            K-M           K-N-M        0
#
# So to compute it, we only need to multiply K-M down to K-M-N, then multiply
# from K down to K-N, then divide the former by the latter. This takes O(N)
# bignum multiplications.

puts "Exact:"
# NOTE: The 1_000_000_000 is a hack to get significant digits after the integer
# division, when neither the numerator or denominator will even *fit* in
# a floating point number. This is probably not the ideal value.
puts 1 - ((mul * 1_000_000_000) / div).to_f / 1_000_000_000.0

puts "1 - (1-G/K)^N estimate:"
puts 1 - (1-g.to_f/k.to_f)**n

puts "GN/K estimate:"
puts g.to_f * n.to_f / k