数概率

【Puzzle】不公平的硬币4

字数统计: 3.4k字   |   阅读时长: 14分
,

摘要: 不公平的硬币4

【对算法，数学，计算机感兴趣的同学，欢迎关注我哈，阅读更多原创文章】
我的网站：潮汐朝夕的生活实验室
我的公众号：算法题刷刷
我的知乎：潮汐朝夕
我的github：FennelDumplings
我的leetcode：FennelDumplings

这是概率面试题连载第 21 期，往期的内容整理在这篇文章里；或者看这个 github 仓库。

目前的题目来源主要是《Fifty challenging problems in probability with solutions》这本书以及朋友在面试中遇到的题。对于每道题，首先给出思路参考，然后用蒙特卡洛模拟的方式验证结果。

在文章 Puzzle-不公平的硬币 中，我们看了一个贝叶斯公式里面比较简单但是很有代表性的一个题。今天我们在那道题的基础上做一点改变形成一个新的题目来看看。为了方便对比，我们首先把上一期的题目抄一遍如下

有 10 枚硬币，其中有 1 枚是不公平的，随机抛出正面的概率为 0.8，另外 9 枚都是公平的，也就是随机抛出正面的概率是 0.5。

现在从 10 枚硬币中随机取出 1 枚，随机地抛了 3 次，结果为”正正反”。求该硬币为不公平硬币的概率。

在文章 Puzzle-不公平的硬币3 中，我们做了一个变种，将【求特定假设成立的概率】改为了【求预测正确的概率】。这个问题稍微有点复杂，在文章 Puzzle-不公平的硬币2 中，我们先处理了一个简化版的问题，并提炼出了这类问题的方法，然后用该方法解决了这个变种问题。

本文我们在以上变种的基础上，再换一种问法，本题由于比较复杂，需要用有向图来对整体流程进行建模。下面我们来看这个题。

问题

有 10 枚硬币，其中有 3 枚是不公平的，随机抛出正面的概率为 0.8，另外 7 枚都是公平的，也就是随机抛出正面的概率是 0.5。

现在从 10 枚硬币中随机取出 1 枚，可以随机地抛 3 次，根据结果来预测该硬币是否为不公平硬币。

目前已经抛出了第1次。

(1) 如果第一次抛出的结果是正面，则三次抛出完成后预测正确的概率是多少；
(2) 如果第一次抛出的结果是反面，则三次抛出完成后预测正确的概率是多少；

(1) 问题对比

本题与文章 Puzzle-不公平的硬币3 解决的问题相比，仅仅是增加了一个已知第一次抛出的结果，求的是已知第一次抛出结果后，预测正确的条件概率。为了方便对比，我们首先把上一期的题目抄一遍如下

有 10 枚硬币，其中有 3 枚是不公平的，随机抛出正面的概率为 0.8，另外 7 枚都是公平的，也就是随机抛出正面的概率是 0.5。

现在从 10 枚硬币中随机取出 1 枚，可以随机地抛 3 次，根据结果来预测该硬币是否为不公平硬币。求预测正确的概率。

该题求的是未知第一次抛出结果时的预测正确的概率，答案是 0.7661，计算过程可以参考上一期文章。

对于本题，我们还是先把假设和证据分别是什么列清楚，然后把各个假设的先验概率以及各个假设下得到各个证据的似然值梳理清楚，再分析已知第一次抛出的结果后，预测正确的条件概率的计算。

(2) 假设和证据

首先我们看都有哪些可能的证据。由于可以抛出 3 次，因此观察到的证据可能有”正正正”、”正正反”、”正反正”、”正反反”, “反正正”、”反正反”、”反反正”、”反反反”这八种。分别记为 e1, e2, e3, e4, e5, e6, e7, e8。

然后我们看都有哪些可能的假设。由于有公平和不公平这两种硬币，因此针对随机取出的硬币有两种假设，”公平” 和 “不公平”，分别记为 H1, H2。

(3) 先验概率

由于本题是 10 枚硬币，其中有 3 枚不公平硬币，因此

(4) 似然值

公平硬币抛出正面的概率为 p = 0.5，不公平硬币抛出正面的概率是 q = 0.8。

如果硬币是公平硬币，则得到 8 种证据的概率均为 1/8，也就是 8 个似然值均为 1/8

如果硬币是不公平硬币，则得到 8 种证据的概率分别为

(5) 第一次抛出后的情况

第一次抛出可能有两种结果：正面和反面，分别记为 X1 和 X2。如果是正面，最终预测对的概率就是 P(T|X1)，如果是反面，最终预测对的概率就是 P(T|X2)。

P(T|X1) 和 P(T|X2) 正是我们要求的值。这两个值在整个获取证据流程中是中间的位置，上游有整体的正确概率 P(T)；下游有三次抛出完成后，各个证据下的正确概率 P(T|e)，还是有点复杂的，我们借助有向图来理解会比较直观。

(6) 有向图建模

X1 对应 e1, e2, e3, e4 这四种证据，因此问当第一次抛出结果为 X1 时预测正确的概率，相当于问 3 次抛出结果后取得的证据为 e1, e2, e3, e4 这四种之一时预测正确的概率。

X2 对应 e5, e6, e7, e8 这四种证据，因此问当第一次抛出结果为 X2 时预测正确的概率，相当于问 3 次抛出结果后取得的证据为 e5, e6, e7, e8 这四种之一时预测正确的概率。

为了理解以上的公式，我们把整个流程抽象成一个有向图，如下：

三层节点的含义如下

  • 根节点(第一层节点)表示一次还没有抛出，还可以抛出 3 次的状态，该状态的值表示整体的预测正确的概率。
  • 第二层节点表示还剩下两次没有抛出的状态，根据第一次的结果可以分为两种状态，状态的值表示在该状态下预测正确的概率。
  • 第三层节点表示三次抛出均完成的状态，根据三次抛出的结果可以分为八种状态(8 种可能证据)，状态的值表示在该证据下下预测正确的概率。

有向边表示从一个状态跳向另一个状态的概率。要计算某节点表示的状态对应的预测正确概率，只需找到该节点的所有出边对应的下一个节点持有的概率值，与边权表示的状态转移概率进行加权求和后就可以得到。

例如在上图中，如果我们要求 P(T)，只需要找到两个出边对应的下游节点 P(T|X1) 和 P(T|X2)，以及转移概率 P(X1) 和 P(X2)，即可得到

如果要求 P(T|X1)，只需要找到四个出边对应的下游节点 P(T|e1), P(T|e2), P(T|e3), P(T|e4)，以及转移概率 P(e1|X1), P(e2|x1), P(e3|X1), P(e4|x1)，即可得到

(7) 推导 P(T|X)

其中 P(ei|Hj) 是我们之前算好的似然值，P(Hj) 是之前算好的先验概率。

由于 X1 包含 ei，即 ei 成立则 X1 一定也成立，所以

因此 P(T|X1) 和Ｐ(T|X2) 的最终计算公式如下

这里 P(X1) 和 P(X2) 需要额外算一下，用到全概率公式

(8) 计算 P(T|X)

将先验概率，似然值，和刚刚算好的 P(X1), P(X2) 代入公式计算

蒙特卡洛模拟

模拟的逻辑与 Puzzle-不公平的硬币3 中的基本一样：模拟的核心是当得到证据 e 之后，如何给出预测结果，也就是代码中 model 实例中的 predict 方法。根据前面的推导，这取决于 16 个似然值。首先随机生成数据模拟出 8 个后验概率(代码中 model 实例的 train 方法)，再随机生成数据模拟预测正确的概率(代码中 model 实例的 test 方法)。只是在预测的时候，增加了对条件概率 P(T|X1) 和 P(T|X2) 的模拟，验证前面的推导过程。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
import numpy as np
from multiprocessing import Pool

class Coin:
    def __init__(self, p, t):
        self.p = p
        self.t = t # 硬币种类 id

    def flip(self):
        """
        如果抛出正面，返回 True
        """
        return np.random.rand() < self.p

class Model:
    def __init__(self, coins, n_flip):
        self.coins = coins
        self.n_flip = n_flip
        # mapping 记录每种证据下两种硬币种类的次数
        # 共有 16 个参数
        self.mapping = {}
        for i in range(1 << n_flip):
            e = [0] * n_flip
            for j in range(n_flip):
                e[n_flip - 1 - j] = i >> j & 1
            self.mapping[tuple(e)] = np.array([0, 0])

    def update(self, e, idx):
        """
        根据数据 (e, idx) 更新 mapping 的参数
        数据的含义是证据 e 对应硬币种类 idx
        """
        self.mapping[e][idx] += 1

    def train(self):
        """
        随机生成数据训练
        """
        n_train = int(1e5)
        for i in range(n_train):
            idx = np.random.randint(0, len(self.coins))
            coin = self.coins[idx]
            e = []
            for _ in range(self.n_flip):
                e.append(int(coin.flip()))
            e = tuple(e)
            self.update(e, coin.t)
        for k, v in self.mapping.items():
            print("证据 {} 下的正常硬币有{}个, 特殊硬币有{}个".format(k, v[0], v[1]))
        print("===训练结束===")

    def predict(self, e):
        """
        根据证据 e 猜硬币种类，返回 0 或 1
        """
        return np.argmax(self.mapping[e])

    def test(self, T):
        correct = 0
        np.random.seed()
        n_X1 = 0
        n_X2 = 0
        n_X1_corrent = 0
        n_X2_corrent = 0
        for _ in range(T):
            # 随机选硬币
            idx = np.random.randint(0, len(self.coins))
            coin = self.coins[idx]
            # 抛出 n_flip 次，获取证据
            e = []
            # 抛出 n_flip 次
            for _ in range(self.n_flip):
                e.append(int(coin.flip()))
            e = tuple(e)
            if e[0] == 1:
                n_X1 += 1
            else:
                n_X2 += 1
            # 根据证据猜硬币种类
            if self.predict(e) == coin.t:
                correct += 1
                if e[0] == 1:
                    n_X1_corrent += 1
                else:
                    n_X2_corrent += 1
        return correct / T, n_X1_corrent / n_X1, n_X2_corrent / n_X2

coins = [Coin(0.5, 0)] * 7 + [Coin(0.8, 1)] * 3
model = Model(coins, 3)
model.train()
T = int(1e7)
args = [T] * 8
pool = Pool(8)
ts = pool.map(model.test, args)
ans = [0, 0, 0]
for t in ts:
    ans[0], ans[1], ans[2] = t
print("P(T): {:.6f}".format(ans[0]))
print("P(T|X1): {:.6f}".format(ans[1]))
print("P(T|X2): {:.6f}".format(ans[2]))

模拟结果

1
2
3
4
5
6
7
8
9
10
11
12
13
证据 (0, 0, 0) 下的正常硬币有8849个, 特殊硬币有216个
证据 (0, 0, 1) 下的正常硬币有8689个, 特殊硬币有908个
证据 (0, 1, 0) 下的正常硬币有8718个, 特殊硬币有988个
证据 (0, 1, 1) 下的正常硬币有8765个, 特殊硬币有3892个
证据 (1, 0, 0) 下的正常硬币有8747个, 特殊硬币有926个
证据 (1, 0, 1) 下的正常硬币有8846个, 特殊硬币有3902个
证据 (1, 1, 0) 下的正常硬币有8709个, 特殊硬币有3810个
证据 (1, 1, 1) 下的正常硬币有8669个, 特殊硬币有15366个
===训练结束===
无额外决策的结果
P(T): 0.766069
P(T|X1): 0.705227
P(T|X2): 0.853590
Share