"""
从目标集合 S 中有放回地抽取 n 次, 计算子集 s 中的每个元素至少抽中 1 次的概率.
目标集合 S 为 [A, A, B, C...] 这种有重复的集合.
集合 s 为 S 的子集, 且 s 中不包含重复元素.
"""
import time
from typing import List
from random import Random
from functools import reduce
def frac(n):
"""
n!
"""
if n <= 1:
return 1
r = 1
while n > 1:
r *= n
n -= 1
return r
def calculate_probability(choices: List[str], targets: List[str], num_draw: int) -> float:
"""
计算此概率:从目标集合 choices 中有放回地抽取 num_draw 次, targets 集合中每个元素至少抽中一个.
targets 必须是 choices 的子集
"""
num_choice = len(choices)
num_target = len(targets)
left_draw = num_draw - num_target
targets = {t: choices.count(t) for t in targets}
if left_draw < 0:
print(f"0.0%")
return 0.0
# 目标集合每个选项出现的次数相乘, 为目标集合的可能数
count = reduce(lambda x, y: x * y, targets.values())
# 如果抽取次数大于目标集合个数, 则还可以再抽取若干次任意元素
count *= (num_choice ** left_draw)
# 以上获得的可能数是固定顺序下的, 所以总的可能数需要再乘以 num_draw 的阶乘
count *= frac(num_draw)
# 每次抽取有 num_choice 次可能, 总共抽取 num_draw 次
total = num_choice ** num_draw
p = round(count / total * 100, 2)
print(f"{p}%")
return p
def observe_probability(choices, targets, num_draw) -> float:
"""
观察此概率:从目标集合 choices 中有放回地抽取 num_draw 次, targets 集合中每个元素至少抽中一个.
targets 必须是 choice 的子集.
"""
rand = Random()
rand.seed(time.time() * 1000)
count = 0
total = 1000000
# 模拟 total 次
for num_test in range(total):
targets = {t: 0 for t in targets}
# 每次模拟抽取 num_draw 次
for draw in range(num_draw):
choice = rand.choice(choices)
if choice in targets:
targets[choice] += 1
values = list(targets.values())
match = all(map(lambda x: x >= 1, values))
if match: # 所有目标都抽到, 计数加 1
count += 1
p = round(count / total * 100, 2)
print(f"{p}%")
return p
def test():
choices = [
"A", "A", "A", "A",
"B", "B", "B", "B",
"C", "C", "C", "C",
"D",
"E",
]
targets = ["A", "B", "C", "D"]
num_draw = 5
observe_probability(choices, targets, num_draw)
calculate_probability(choices, targets, num_draw)
if __name__ == "__main__":
test()
1
ConradG 2021-01-28 22:08:19 +08:00
推测是结果偏大?
如果是,问题出在“如果抽取次数大于目标集合个数, 则还可以再抽取若干次任意元素”。 如果把四个 A 依次标记为 a1 到 a4,则 ABCD+A 这种形式按照 LZ 的计算方法会出现 a1BCD+a2 和 a2BCD+a1 的重复 |
2
necomancer 2021-01-28 22:30:48 +08:00
有标准答案吗?假定 s 是 s 集合的长度,且 P_{i=1,2,...,s} 是抽到 s 中每个元素的概率,那么当 n < s 的时候概率是 0,当 n > s 的时候是 P_{n}^s \Pi_{i=1}^s p_i
|
3
necomancer 2021-01-28 22:31:07 +08:00
n>=s
|
4
vance123 2021-01-28 22:37:24 +08:00 via iPad
不就是随机播放歌单,平均多少次才能听完所有歌的变体吗,我记得期望好像是 7 点几
|
5
necomancer 2021-01-28 22:50:04 +08:00
我感觉模拟可以这么做,反正是有放回抽样,那每次抽到某个元素的概率都是给定的,所以直接当 n>=s 的时候就从一个 s+1 的数组 S=s+np.inf (随便一个什么东西)按概率 ps={p_1,p_2,....,p_s,1-\sum_i^s p_i},最后一个元素只要和其他的都不一样就行了,用
out = np.random.choice(S, size=n,p=ps),然后用个 histogram 或者 hash 判定一下? |
6
gulu OP |
7
gulu OP @necomancer
我感觉模拟的结果是正确的,应该是计算的方法出错了吧 |
8
ConradG 2021-01-28 23:19:43 +08:00 1
@gulu
因为你在最后一步“count *= frac(num_draw)”算了一次全排列,那么就要求之前的统计结果是无序的。而我提到的重复是出现在第二步结束时候。 顺便这类问题的常规解法是写成递归。 |
9
gulu OP |
10
necomancer 2021-01-29 01:48:12 +08:00 1
我好像也算大了.... 等概率的话还有办法的,可以表达成 s! S_2(n, s),我记得这个叫优惠券问题。
|
11
necomancer 2021-01-29 02:53:45 +08:00
p(n,s) = \sum_{j=0}^{s-1} (-1)^{s-1-j}C_{m-j-1}^{m-c}\sum_{|J|=j}P_J^{n-1}(1-P_J)....(1)
其中 P_J 是子集概率加和,P_J := \sum_{i\in J} p_i, \sum_{|J|=j} 代表对所有 |J|=j 大小的子集的加和,m 代表总独立元素个数。对均匀分布来说, p(n,s) = m! S_2(n-1, s-1)/((m-s)! m^n)....(2) 这么看好像模拟的时候不能把 m 个独立元素中选择 s 个折合成 s+1 里选 s 个…… |
12
necomancer 2021-01-29 02:54:23 +08:00
Google 了一下 Coupon collector's problem 发现这玩意儿好像还很难
|
13
necomancer 2021-01-29 02:55:18 +08:00
基础知识的话得看看停时和鞅。
|
14
gulu OP @necomancer
谢谢大佬,我去学习一下 |
15
gulu OP @necomancer
老哥,我初步看了一下“集券”问题,感觉和我这个问题有点不一样。 集券问题: Given n coupons, how many coupons do you expect you need to draw with replacement before having drawn each coupon at least once? 有 n 张优惠券,有放回地抽取多少次,才能保证每张券抽取到至少一次? 我的问题: 有 N 张优惠券,不平均地分类为 n 种。有放回地抽取多少次,才能保证 m 种优惠券至少抽取到一次? |
16
HelloViper 2021-01-29 10:09:06 +08:00
放回取 n 次 == len(S)^n ,
直接 product 出笛卡尔积,然后遍历 set 后 s 是否 in 就结了 |
17
necomancer 2021-01-29 14:27:16 +08:00 1
你 Google 一下 coupon collector's problem with unequal probability of subset?
@gulu |
18
gulu OP @necomancer
这就去看看😁 |