【Puzzle】用户一次游戏带来的收入

字数统计: 4.1k字 | 阅读时长: 17分

2021-05-28

摘要: 来自朋友的面试的概率题

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

写在前面

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

今天依然是看一道朋友面试时遇到的题，我认为本题是概率计算问题中最有代表性的一道题，可以延伸出很多变种，所以准备多交代一些相关内容，以后遇到这类题，90% 都可以不用怕了。

实际上在前面的概率面试题连载文章中，我们解决过类似的问题，只是那道题很简单，因此直接提了思路就直接过了，没有做方法论的总结。

本文会先以此前解决过的这个简单问题为例子，交代一些相关的背景，提炼出一类题的方法论，内容涉及到全期望公式，有向图，期望DP，高斯消元。

然后再引入今天这道题并解决，然后再看一个在今天的题目基础上的变种问题。

$1 回顾: 试验直到第一次成功

回顾在之前的文章中解决的一个问题，问题和解法如下，原文可以看【Puzzle】试验直到第一次成功

一枚骰子掷出第一个 6 平均需要掷多少次。

在那篇文章中我们用两种方法解决了这个问题:

期望的定义

第一种是直接用期望的定义求期望

$E(X) = \sum\limits_{i}i \times p(i)$

由于 i 的所有可能取值和 p(i) 我们都可以算出来，于是通过定义再加上一些公式推导技巧我们求出了 E(X) 的解析解。

全期望公式

第二种是基于【全期望公式】求期望

$\begin{align} E_{X}(X) &= E_{Y}(E_{X|Y}(X|Y)) \\ &= \sum\limits_{y}p(Y=y) \times E_{X|Y=y}(X|Y=y) \\ \end{align}$

当时在那篇文章中的思路是这么写的

记 m 为第一个 6 需要的次数
当第一次投掷成功时，则需要的平均次数是 1
当第一次投掷失败时，则需要的平均次数是 1 + m
于是有以下公式，进而可以求出 m
m = p * 1 + q * (1 + m)

从全期望公式的角度看

X 是掷出第一个 6 的次数，m 就是 E(X)；
Y 是第一次投掷出的点数是否为 6，是则记 Y=1，否则记 Y=0；
Y=1 的概率为 p，Y=0 的概率为 q；
E(X|Y=1) = 1，E(X|Y=0) = 1 + E(X)

代入公式即可得到解题思路中用的公式 m = p * 1 + q * (1 + m)

$\begin{align} E(X) &= p(Y=1)E(X|Y=1) + p(Y=0)E(X|Y=0) \\ &= p \times 1 + q \times (1 + E(X)) \\ \end{align}$

有向图

用全期望公式的求期望的过程可以用有向图的方式表示出来，例如对于上面的问题，我们可以画出下面的有向图

其中节点 0 表示未投掷出 6 的状态，节点 1 为表示掷出 6 的状态。0 为起点，1 为终点。

图中的每条有向边有一个概率值，还有另一个权值，节点有一个期望值。

我们要求的是从起始状态走向终点状态要走的步数的期望，每个节点有一个期望值表示从该节点走到终点节点的期望。由于 1 是终点，所以 E[1] = 0。0 是起点，我们要求的就是 E[0]。

从节点 0 会伸出两条边，一条回到 0，表示当前的投掷没有掷出 6，另一条通向 1，表示当前的投掷掷出了 6。边上的概率值表示当前投掷的各个值的概率，而投掷一次步数加一，因此额外的边权都是 1。

定义好这个图之后，我们就可以写出 E[0] 的表达式了

1	E[0] = p * (1 + E[1]) + q * (1 + E[0])

这个表达式就是全期望公式，只是将全期望公式中 E(X|y) 的部分用图中点权和边权的方式定义了计算方法。

期望DP

通过以上的分析可以发现，当根据具体问题定义好有向图 D，起点 s，终点 t 后，我们就可以用类似于动态规划的方式求从 s 到 t 的某个指标的期望。

状态定义
dp[u] := 从 u 到 t 的某个指标的期望

初始化
dp[t]

答案
dp[s]

状态转移
dp[u] = g(node[u], sum(p[u][v] * f(w[u][v], dp[v])))

其中 f 是边权 w[u][v] (在上面的题中是 1) 和下一点的期望 dp[v] 的函数，g 是当前点的点权 node[u] (在上面的题中没有) 和 sum(p[u][v] * f) 的函数, 具体的函数形式取决于问题(在上面的题中 f 就是简单的 w[u][v] + dp[v], g 没有)，需要具体分析。

求解 dp[s] 的过程：对反图进行拓扑排序，按照拓扑排序的顺序依次求解各个节点 u 的 dp[u]，直至求解到 dp[s]

高斯消元

当建出的有向图中有环的时候，求解 E[s] 的过程如果直接用期望 DP 从 dp[t] 向 dp[s] 推导的话是不行的，因为在推导 DP 状态时会出现类似于下面的情况(就是后面的【用户一次游戏带来的收入】那道题的转移方程，这里提前看一下)

dp[0] = dp[1]
dp[1] = 1 + 0.7 * dp[2] + 0.3 * dp[0]
dp[2] = 2 + 0.6 * dp[3] + 0.4 * dp[1]
dp[3] = 3 + 0.5 * dp[4] + 0.5 * dp[2]
dp[4] = 0

这种情况 dp 状态的转移形成了环，比如要求 dp[1] 要先知道 dp[2]，要求 dp[2] 就要先知道 dp[1]，没法求了。

如果方程组是线性方程组，还有有办法的，解决办法是利用高斯消元的方式整体解这个线性方程组。

关于高斯消元的推导、代码、题目，可以参考这篇文章，高斯消元与线性方程组

$2 方法总结

通过前面的分析，我大致可以总结出求解这类求期望的问题的方法论。

首先理论基础就是全期望公式。

然后针对具体的问题，我们先分析出建图需要的信息

起点 s 和终点 t 状态是什么，以及还有哪些状态节点
求出各个状态节点走一步可以走到哪些状态节点，概率分别多少
分析计算所求的期望是否需要额外的点权和边权

然后建图，初始化 dp[t]

分析状态转移方程中的 g, f 并写出实际的转移方程

最后求解 dp[s]:

如果图没有环，直接按拓扑排序的顺序进行动态规划的方式推导即可；
如果图有环，就看方程组是什么样的，线性方程组可以用高斯消元法求解，其它形式的方程我暂时也不清楚解法。

$3 题目: 装备升级

下面开始今天要看的朋友在面试时遇到的题目。题目描述如下:

玩家的装备有从 0 到 4 这 5 个等级，每次升级都需要一个工序，一共有 0->1, 1->2, 2->3, 3->4 这 4 个工序，成功率分别为 0.9, 0.7, 0.6, 0.5。工序成功，玩家的装备就会升一级，但如果失败，玩家的装备就要倒退一级。例如玩家当前的等级为 2，目前执行 2->3 这个工序，如果成功，则装备等级从 2 升为 3，如果失败，装备等级就从 2 降到 1。

问：玩家装备初始等级为 0, 升到 5 平均要经历多少个工序。

思路参考

建图

按照前面的方法总结。我们首先分析题目并建图，点状态是 0，终点状态是 4，此外还有 1, 2, 3 这三个中间状态。状态间的转移方向以及概率题目中明确给了。我们要求的是从 s=0 走到 t=4 平均需要走的步数。

由于每经过一条边都相当于走了一步，所以边上有额外的权 1。除此之外没有别的权，节点上也没有权。根据这些信息，我们首先把图建出来，如下

期望DP

期望 DP 的方程组如下:

dp[0] = 0.9 * (1 + dp[1]) + 0.1 * (1 + dp[0])
dp[1] = 0.7 * (1 + dp[2]) + 0.3 * (1 + dp[0])
dp[2] = 0.6 * (1 + dp[3]) + 0.4 * (1 + dp[1])
dp[3] = 0.5 * (1 + dp[4]) + 0.5 * (1 + dp[2])
dp[4] = 0

这是一个有环的图，且方程是线性方程组，面试的时候手推就可以。

如果要编程的话，需要整理成标准形式后用高斯消元解决。

高斯消元

标准形式如下

$\begin{gathered} \begin{bmatrix} 0.9 & -0.9 & 0 & 0 & 0 \\ -0.3 & 1 & -0.7 & 0 & 0 \\ 0 & -0.4 & 1 & -0.6 & 0 \\ 0 & 0 & -0.5 & 1 & -0.5 \\ 0 & 0 & 0 & 0 & 1 \\ \end{bmatrix} \begin{bmatrix} dp[0] \\ dp[1] \\ dp[2] \\ dp[3] \\ dp[4] \\ \end{bmatrix} = \begin{bmatrix} 1 \\ 1 \\ 1 \\ 1 \\ 0 \\ \end{bmatrix} \end{gathered}$

高斯消元求解

#include <vector>
#include <cmath>
#include <iostream>
#include <iomanip>

using namespace std;

int main()
{
    // 系数矩阵
    vector<vector<double>> c{{0.9, -0.9, 0, 0, 0}
                            ,{-0.3, 1, -0.7, 0, 0}
                            ,{0, -0.4, 1, -0.6, 0}
                            ,{0, 0, -0.5, 1, -0.5}
                            ,{0, 0, 0, 0, 1}
                            };
    // 常数列
    vector<double> b{1, 1, 1, 1, 0};
    int n = 5;
    const double EPS = 1e-10;

    // 高斯消元, 保证有唯一解
    for(int i = 0; i < n; ++i)
    {
        // 找到 x[i] 的系数不为零的一个方程
        for(int j = i; j < n ;++j)
        {
            if(fabs(c[j][i]) > EPS)
            {
                for(int k = 0; k < n; ++k)
                    swap(c[i][k], c[j][k]);
                swap(b[i], b[j]);
                break;
            }
        }
        // 消去其它方程的 x[i] 的系数
        for(int j = 0; j < n; ++j)
        {
            if(i == j) continue;
            double rate = c[j][i] / c[i][i];
            for(int k = i; k < n; ++k)
                c[j][k] -= c[i][k] * rate;
            b[j] -= b[i] * rate;
        }
    }
    cout << std::fixed << std::setprecision(6);
    for(int i = 0; i < n; ++i)
    {
        cout << "dp[" << i << "] = " << b[i] / c[i][i] << endl;
    }
}

求解结果

dp[0] = 10.888889
dp[1] = 9.777778
dp[2] = 7.873016
dp[3] = 4.936508
dp[4] = 0.000000

蒙特卡洛模拟验证结果

import numpy as np
import bisect

class Game:
    def __init__(self, p):
        self.n = p.shape[0]
        self.p = p
        for i in range(self.n):
            for j in range(1, self.p.shape[1]):
                self.p[i][j] += self.p[i][j - 1]

    def step(self, pos):
        return bisect.bisect_left(self.p[pos]
                                 ,np.random.rand())

    def __call__(self):
        pos = 0
        n_steps = 0
        while pos != self.n:
            pos = self.step(pos)
            n_steps += 1
        return n_steps

p = np.array([[0.1, 0.9, 0.0, 0.0, 0.0]
             ,[0.3, 0.0, 0.7, 0.0, 0.0]
             ,[0.0, 0.4, 0.0, 0.6, 0.0]
             ,[0.0, 0.0, 0.5, 0.0, 0.5]
             ])
game = Game(p)

def test():
    T = int(1e7)
    total_n_steps = 0
    for t in range(T):
        total_n_steps += game()
    print("Average Steps: {:.6f}".format(total_n_steps / T))

for i in range(10):
    test()

模拟结果

Average Steps: 10.890664
Average Steps: 10.887758
Average Steps: 10.882893
Average Steps: 10.889820
Average Steps: 10.891342
Average Steps: 10.888397
Average Steps: 10.888496
Average Steps: 10.891009
Average Steps: 10.890464
Average Steps: 10.886080

$4 变种问题: 一次游戏通关带来的收入

下面我们在装备升级问题上做一些变化，我们保持图结构不变，但是我们把期望 DP 转移方程中的边权改为节点的权。

我们考虑一个在游戏相关岗位面试常见的问题: 一次游戏通关带来的收入

一个游戏有四关，通过概率依次为0.9, 0.7, 0.6, 0.5。
第一关不收费，第二到四关每次收费分别为1块, 2块, 3块。
用户每玩一次都会无限打下去直至通关，通关后用户可以提现 10 块钱作为奖励。
问: 公司可以在每次用户游戏中平均挣多少钱。

思路参考

我们首先考虑【公司可以在每次用户游戏中收费多少钱】，然后减去 10 块钱的奖励就是挣的钱。

建图

图与上面的装备升级那道题一样，只是边权没有了，节点有权重表示费用。

期望 DP

期望 DP 的方程组如下:

dp[0] = dp[1]
dp[1] = 1 + 0.7 * dp[2] + 0.3 * dp[0]
dp[2] = 2 + 0.6 * dp[3] + 0.4 * dp[1]
dp[3] = 3 + 0.5 * dp[4] + 0.5 * dp[2]
dp[4] = 0

这是一个有环的图，且方程是线性方程组，面试的时候手推就可以。

如果要编程的话，需要整理成标准形式后用高斯消元解决。

高斯消元

标准形式如下

$\begin{gathered} \begin{bmatrix} 1 & -1 & 0 & 0 & 0 \\ -0.3 & 1 & -0.7 & 0 & 0 \\ 0 & -0.4 & 1 & -0.6 & 0 \\ 0 & 0 & -0.5 & 1 & -0.5 \\ 0 & 0 & 0 & 0 & 1 \\ \end{bmatrix} \begin{bmatrix} dp[0] \\ dp[1] \\ dp[2] \\ dp[3] \\ dp[4] \\ \end{bmatrix} = \begin{bmatrix} 0 \\ 1 \\ 2 \\ 3 \\ 0 \\ \end{bmatrix} \end{gathered}$

高斯消元求解

#include <vector>
#include <cmath>
#include <iostream>
#include <iomanip>

using namespace std;

int main()
{
    // 系数矩阵
    vector<vector<double>> c{{1, -1, 0, 0, 0}
                            ,{-0.3, 1, -0.7, 0, 0}
                            ,{0, -0.4, 1, -0.6, 0}
                            ,{0, 0, -0.5, 1, -0.5}
                            ,{0, 0, 0, 0, 1}
                            };
    // 常数列
    vector<double> b{0, 1, 2, 3, 0};
    int n = 5;
    const double EPS = 1e-10;

    // 高斯消元, 保证有唯一解
    for(int i = 0; i < n; ++i)
    {
        // 找到 x[i] 的系数不为零的一个方程
        for(int j = i; j < n ;++j)
        {
            if(fabs(c[j][i]) > EPS)
            {
                for(int k = 0; k < n; ++k)
                    swap(c[i][k], c[j][k]);
                swap(b[i], b[j]);
                break;
            }
        }
        // 消去其它方程的 x[i] 的系数
        for(int j = 0; j < n; ++j)
        {
            if(i == j) continue;
            double rate = c[j][i] / c[i][i];
            for(int k = i; k < n; ++k)
                c[j][k] -= c[i][k] * rate;
            b[j] -= b[i] * rate;
        }
    }
    cout << std::fixed << std::setprecision(6);
    for(int i = 0; i < n; ++i)
    {
        cout << "dp[" << i << "] = " << b[i] / c[i][i] << endl;
    }
}

求解结果

dp[0] = 16.000000
dp[1] = 16.000000
dp[2] = 14.571429
dp[3] = 10.285714
dp[4] = 0.000000

因此公司可以在每次用户游戏中收费 dp[0] = 16 块钱，减去用户通关的 10 块钱奖金，公司可以挣 6 块钱。

蒙特卡洛模拟验证结果

import numpy as np
import bisect

class Game:
    def __init__(self, transfer, cost):
        self.n = transfer.shape[0]
        self.cost = cost
        self.transfer = transfer
        for i in range(self.n):
            for j in range(1, self.transfer.shape[1]):
                self.transfer[i][j] += self.transfer[i][j - 1]

    def step(self, pos):
        return bisect.bisect_left(self.transfer[pos]
                                 ,np.random.rand())

    def __call__(self):
        pos = 0
        pay = 0
        while pos != self.n:
            pay += self.cost[pos]
            pos = self.step(pos)
        pay += self.cost[self.n]
        return pay

transfer = np.array([[0.1, 0.9, 0.0, 0.0, 0.0]
                    ,[0.3, 0.0, 0.7, 0.0, 0.0]
                    ,[0.0, 0.4, 0.0, 0.6, 0.0]
                    ,[0.0, 0.0, 0.5, 0.0, 0.5]
                    ])
cost = np.array([0.0, 1.0, 2.0, 3.0, -10.0])
game = Game(transfer, cost)

def test():
    T = int(1e5)
    total_pay = 0
    for t in range(T):
        total_pay += game()
    print("Average income: {:.4f}".format(total_pay / T))

for i in range(10):
    test()

模拟结果

Average income: 5.9768
Average income: 6.0043
Average income: 5.9052
Average income: 6.0664
Average income: 6.0720
Average income: 5.9687
Average income: 6.0101
Average income: 6.0450
Average income: 5.9871
Average income: 6.0704

潮汐朝夕的生活实验室 \Doge 陪伴一个算法工程师的职业生涯

数概率

【Puzzle】用户一次游戏带来的收入

写在前面

$1 回顾: 试验直到第一次成功

期望的定义

全期望公式

有向图

期望DP

高斯消元

$2 方法总结

$3 题目: 装备升级

思路参考

建图

期望DP

高斯消元

蒙特卡洛模拟验证结果

$4 变种问题: 一次游戏通关带来的收入

思路参考

建图

期望 DP

高斯消元

蒙特卡洛模拟验证结果