三国杀武将技能概率
前言
本文的主体内容取自于本人十年前在百度帖吧三国杀吧的一篇帖子：【理论】数学与三国杀——三国杀武将技能概率。现加以整理、改进和删增。
本文通过研究《三国杀》中部分武将技能所涉及的概率与期望，对这些技能的强度与收益进行定量刻画，从而为玩家在特定对局中的选将、技能发动时机及出牌顺序提供决策参考。
研究范围
卡牌
本文只考察军争篇卡牌。其点数、花色以及类型具体分布如下：
（图片来源：https://www.cnblogs.com/yhm138/articles/14091928.html，有勘误）
（图片来源：https://www.cnblogs.com/yhm138/articles/14091928.html）
数学模型
为便于统一讨论，本文将技能效果抽象为随机变量，并对其概率分布及数学期望进行分析。文中所称“期望”，均指技能触发结果或收益函数所对应随机变量的数学期望。
理想模型
若无特殊说明，本文所采用的数学模型均为理想模型，即假设任意一个时刻，牌堆中的卡牌数量无限多，且某牌按照其在上述军争篇中的概率均匀分布（并假设各次抽牌在概率意义上相互独立）。例如：
- 任一花色的概率为$1/4$
- 任一点数的概率为$1/13$1
- 普通（非延时类）锦囊牌的概率为$43/160$
该假设并非对真实牌堆的物理描述，而是用于简化计算、刻画长期平均行为的数学近似。
无限牌堆模型实际上意味着有放回抽样。
在此基础上，关于拼点，我们进一步假设拼点双方亮出的牌的点数均匀分布。由对称性或直接计数易知，拼点赢的概率为$\frac{6}{13}$。
此外，本文中所有的小数结果仅保留四位有效数字。
现实模型
在某些情况下，本文也会采用现实模型。该模型有以下特点：
- 游戏模式为军争篇八人场
- 不放回抽样
- 假设在整局游戏过程中：
- 场上平均存活角色数为 4
- 每名角色平均拥有手牌 2 张、装备牌 2 张
- 场上平均存在 1 张判定牌。
这样，牌堆中的牌数平均为：$160-4*(2+2)-1=143$。
理想模型用于刻画技能的统计平均强度与内在概率结构，而现实模型则用于评估这些结论在标准军争局中的数值偏差。
武将
文中所涉及到的武将选取自以下系列：
- 标准版
- 神话再临·风火林山阴雷
- 一将成名2011-2015
- 原创之魂2016-2017
- 一将之魂
- SP
另外，有些武将技能的概率计算过于简单直白，也不在本篇的讨论之列，如:
- 夏侯惇刚烈（成功概率$3/4$）
- 神关羽武魂（失败概率$13/160 = 0.08125$）
- 顾雍秉壹（摸牌概率$1/2^{X-1}$）
- 辛宪英（原创之魂2017）忠鉴（手牌上限减少概率$1/2*12/13 = 6/13 \approxeq 0.4615$）。
- 张星彩枪舞（参考下表【杀】的点数分布）
| 点数 | A | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | J | Q | K |
| 【杀】数量 | 0 | 1 | 1 | 4 | 4 | 4 | 6 | 7 | 5 | 8 | 3 | 0 | 1 |
- 秦宓天辩（拼点赢概率可拆解为若干互斥情形的加权和，其计算结果为：$3/52*1/13\sum_{i=0}^{11}i+16/52 = 203/338 \approxeq 0.6006$）
上述技能要么概率结构过于直接，要么计算不涉及复杂的条件依赖，因此不作为本文重点讨论对象。
武将技能概率
伯努利试验
本文采用负二项分布“在达到固定成功次数\(r\)前所经历的失败次数（或总试验次数）”的参数化方式。
有一些武将技能概率问题符合典型的伯努利试验下的概率分布。我们以邓艾、甄姬和孟获三名武将的技能为例一一说明。对每名武将都会采用理想模型（有放回）和现实模型（不放回）进行分析。
邓艾（神话再临·山）
屯田：每次当你于回合外失去牌时，可进行一次判定，将非♥结果的判定牌置于你的武将牌上，称为“田”；每有一张田，你计算与其他角色的距离便-1。
凿险：觉醒技，回合开始阶段，若田的数量达到3张或更多，你须减1点体力上限，并永久获得技能“急袭”（出牌阶段，你可以把任意一张田当【顺手牵羊】使用）。
下面我们来研究可以发动觉醒技凿险时所需的屯田次数。
【理想模型】
定义发动觉醒技凿险时所需的屯田次数为$X$，并定义一次伯努利试验（Bernoulli trial）的“成功”为判定非♥，而“失败”为判定♥。由此可知，$X$为成功和失败次数之和。若定义成功次数为$r$，则$X$符合负二项分布（Negative binomial distribution）。
概率质量函数（PMF, Probability mass function）:
\[p_{X}(n) = {\binom{n-1}{k}}p^{r}(1-p)^{k}, \quad n = k+r \enspace \text{and} \enspace k \in \{\,0,1,2,\dots\,\}\]注：此处采用负二项分布的“以失败概率$p$入式”的参数化；即$p$表示判定为♥的概率，$r$表示达到$r$次非♥所需总试验次数。
在此模型下，$r=3, \enspace p=3/4$。
其概率曲线为：
重点观察$n$较小时的概率散点图：
数学期望（Expected value）和方差（Variance）分别为：
\[\boxed{ \mathbb E[X] = \frac{r}{p} = 4 }\] \[\boxed{ \operatorname{Var}(X) = \frac{r(1-p)}{p^{2}} = \frac{4}{3} }\]这说明在理想情况下，发动觉醒技凿险时所需的屯田次数平均为4次。
【现实模型】
定义发动觉醒技凿险时所需的屯田次数为$X$，牌堆中的牌数为$N$，牌堆中♥的牌数为$K$ ($K \le N$)，并定义一次伯努利试验的“成功”为判定♥，而“失败”为判定非♥（注意：此处为方便使用负超几何分布的标准形式，成功与失败的定义与理想模型下相反）。若定义失败次数为$r$，则成功次数符合负超几何分布（Negative hypergeometric distribution）。
PMF：
\[p_{X}(k)=\frac{\binom{K+r-1}{k}\binom{N-r-K}{K-k}}{\binom{N}{K}}, \quad k=0,1,\dots,K\]在此模型下，$r=3$。
数学期望和方差分别为：
\[\boxed{ \mathbb E[X] = \frac{rK}{N-K+1}+r = \frac{rN+r}{N-K+1} = \frac{3N+3}{N-K+1} }\] \[\boxed{ \operatorname{Var}(X) = \frac{rK(N+1)(N-K-r+1)}{(N-K+1)^{2}(N-K+2)} = \frac{3K(N+1)(N-K-2)}{(N-K+1)^{2}(N-K+2)} }\]如果我们假定$K = N/4$，那么：
\[\mathbb E[X] = 4-\frac{4}{3N+4} < 4\] \[\operatorname{Var}(X) = \frac{12N(N+1)(3N-8)}{(3N+4)^{2}(3N+8)}\]可见，$\mathbb E[X]$与$N$正相关，且小于理想模型下屯田次数的数学期望。
甄姬（标准版）
洛神：回合开始阶段，你可以进行判定：若为♠或♣，立即获得此牌，并可以再次使用洛神——如此反复，直到出现♥或♦为止。
下面我们考虑洛神所获得的牌数。
【理想模型】
定义洛神获得的牌数为$Y$，并定义一次伯努利试的“成功”为判定红牌，而“失败”为判定黑牌。很显然，$Y$属于失败次数，且其概率分布符合几何分布（Geometric distribution）的第二种形式。
PMF：
\[p_{Y}(k) = (1-p)^{k}p, \quad k \in \{\,0,1,2,\dots\,\}\]其中$p$为“失败”的概率，在此模型下$p = 1/2$。
其概率曲线为：
重点观察$k$较小时的概率散点图：
数学期望和方差分别为：
\[\boxed{ \mathbb E[Y] = {\frac{1-p}{p}} = 1 }\] \[\boxed{ \operatorname{Var}(Y) ={\frac {1-p}{p^{2}}} = 2 }\]这说明在理想情况下，洛神所得的牌数平均为1张。
实际上，几何分布是负二项分布的特殊情况：
\[\operatorname{Geom}(p) = \operatorname{NB}(1,p)\]我们由此可以发现甄姬洛神和邓艾屯田的内在联系。
接下来我们考察甄姬的Online战功“洛神赋”：使用甄姬一回合内发动洛神在不被改变判定牌的情况下连续判定黑色花色至少8次。
\[p_{Y}(k \ge 8) = 1 - \sum_{k = 0}^7{p_Y(k)} = \frac{1}{256} \approxeq 0.003906\]【现实模型】
定义洛神获得的牌数仍为$Y$，牌堆中的牌数为$N$，牌堆中红牌数为$K$ ($K \le N$)，并定义一次伯努利试验的“成功”为判定黑牌，而“失败”为判定红牌（注意：此处的定义与理想模型下刚好相反）。我们已经知道几何分布是负二项分布的特例，所以成功次数符合负超几何分布（Negative hypergeometric distribution），且失败次数$r=1$。
数学期望和方差分别为：
\[\boxed{ \mathbb E[Y] = \frac{rK}{N-K+1} = \frac{K}{N-K+1} }\] \[\boxed{ \operatorname{Var}(Y) = \frac{rK(N+1)(N-K-r+1)}{(N-K+1)^{2}(N-K+2)} = \frac{K(N+1)(N-K)}{(N-K+1)^{2}(N-K+2)} }\]我们也可以通过Stars and bars的方法来直观地得到上述数学期望。设想有$N$个球，其中$K$个为黑色，那么红色球有$N-K$个。把这些红球看作隔板，它们把空间分成了$N-K+1$个区域。每个黑球在某一个区域的概率都相同，为$\frac{1}{N-K+1}$。因此，$K$个黑球在第一个区域的期望为$\frac{K}{N-K+1}$。
如果我们假定$K = N/2$，那么：
\[\mathbb E[Y] = \frac{1}{1+\frac{2}{N}} < 1\] \[\operatorname{Var}(Y) = \frac{N^2(N+1)}{(N+2)(N+4)}\]可见，$\mathbb E[Y]$与$N$正相关，且小于理想模型下洛神牌的数学期望。
前文对军争篇模型的介绍中已经提到，我们假设某一时刻，牌堆平均有143张牌。那么：
\[\mathbb E_{\text{avg}} = \frac{1}{1+\frac{2}{143}} \approxeq 0.9862\]孟获（神话再临·林）
再起：摸牌阶段，若你已受伤，你可以放弃摸牌并展示牌堆顶的X张牌，X为你已损失的体力值，其中每有一张♥牌，你回复1点体力，然后弃掉这些♥牌，将其余的牌收入手牌。
接下来我们来研究再起获得的体力值。
【理想模型】
定义在损失$X \enspace (0 < X \le 4)$体力值的情况下，再起所获得的体力值为$H$，并定义一次伯努利试验的“成功”为判定♥，而“失败”为判定非♥。很显然，$H$属于成功次数，且其概率分布符合二项分布（Binomial distribution）$H \sim \operatorname{Bin}(n,p)$。
PMF：
\[p_{H}(k) = {\binom{n}{k}}p^{k}(1-p)^{n-k}, \quad k \in \mathbb{N} \enspace \text{and} \enspace k=0,1,2,\dots,n\]其中$n$为伯努利试验的次数，$p$为判定为♥的概率。在此模型下，$n=X, \enspace p=1/4$。
其概率散点图为：
以$X = 3$为例：
数学期望和方差分别为：
\[\boxed{\mathbb E[H] = \frac{X}{4}}\] \[\boxed{\operatorname{Var}(H) = np(1-p) = \frac{3X}{16}}\]同理，定义所获得手牌为$C$，则$C \sim \operatorname{Bin}(X,3/4)$。
\[\boxed{\mathbb E[C] = \frac{3X}{4}}\] \[\boxed{\operatorname{Var}(C) = np(1-p) = \frac{3X}{16}}\]这里我们继续讨论一个问题：$X=?$时发动再起收益最高？
我们需要首先定义“收益”的概念。如果认为一点体力值等价于$R$张牌，那么我们可以将再起的收益用牌数来表示：
\[P = C + HR\] \[\mathbb E[P] = \mathbb E[C]+R\mathbb E[H] = \frac{X}{4}(3+R)\]一般认为，$R \approxeq 2$，即“一血=两牌”理论（郭嘉遗计、黄盖苦肉、钟会排异等都是佐证）。该理论经大量游戏实践验证而一度在玩家中非常流行。当然，不少玩家对此理论也持不同观点，因此我们要讨论在再起技能的情况下，$R$值的合理范围。
如果不发动再起，摸牌阶段会获得$2$张牌。若要使发动再起的收益不劣于不发动，那么：
\[\mathbb E[P] \ge 2\]所以
\[R \ge \frac{8}{X} - 3\]上述曲线图如下：
可见：
- 当$X \ge 3$时，无论$R$值如何，发动再起都有正收益
- 当$X = 1$时，我们应当认为发动再起应当有负收益，否则玩家可以在该角色受伤时总是选择发动此技能而不会担心带来负收益，而这就意味着此技能设计过强。所以，$R < 5$
- 当$X = 2$时，如果我们认为发动再起应当有正向收益，那么$R \ge 1$；否则$0 < R \le 2$
综上所述，$R = 2$是一个合理的取值。在此取值下，当$X le 2$时，若非特殊情况，玩家不应发动此技能。
【现实模型】
定义牌堆中的牌数为$N$，牌堆中♥牌数为$K$ ($K \le N$)。在此模型下，伯努利试验是不放回的取样，此时成功次数符合超几何分布（Hypergeometric distribution）$H \sim HG_{N,K,n}(k)$。
PMF：
\[p_{H}(k)=\frac{\binom{K}{k}\binom{N-K}{n-k}}{\binom{N}{n}}, \quad k \in \{\,0,1,2,\dots\,\} \enspace \text{and} \enspace \max(0,n+K-N)\leq k\leq \min(K,n)\]其中$n$为伯努利试验的次数，$n=X$。
数学期望和方差分别为：
\[\boxed{\mathbb E[H] = np = Xp}\] \[\boxed{\operatorname{Var}(H) = \frac{np(1-p)(N-n)}{N-1} = \frac{Xp(1-p)(N-X)}{N-1}}\]其中，$p = K/N$。
如果我们假定$K = N/4$，则$p=1/4$，那么：
\[\mathbb E[H] = \frac{X}{4}\] \[\operatorname{Var}(H) = \frac{3X}{16}\frac{N-X}{N-1} < \frac{3X}{16}\]实际上，超几何分布和负超几何分布有以下联系：
\[NHG_{N,K,r}(k)=1-HG_{N,N-K,k+r}(r-1)\]我们由此可以发现孟获再起和邓艾屯田的内在联系。
伯努利试验概率分布小结
以上武将技能分析涵盖了抽样问题的四种概率分布：
| 有放回（理想模型） | 不放回（军争模型） | |
| 抽样数固定时成功数的分布 | 二项分布 再起 | 超几何分布 再起 |
| 失败数固定时成功数的分布 | 负二项分布 屯田（特例：几何分布 洛神） | 负超几何分布 屯田 洛神 |
还有几名武将的技能符合以上几种分布，我们简略列举如下。由于分析过程类似，我们将只讨论理想模型而忽略现实模型。
黄月英（标准版）
集智：当你使用普通锦囊牌时，你可以摸一张牌。
假设当前手牌中的普通锦囊牌有$r$张，并定义本回合内可发动集智的次数为$X$（假设每当使用普通锦囊牌时都会发动集智；不考虑【无中生有】和【五谷丰登】的影响），则$X$符合负二项分布：
\[X \sim NB(r,p)\]其中$p$为判定为非普通锦囊牌的概率。在此模型下，$p=\frac{117}{160}$。
数学期望和方差分别为：
\[\boxed{ \mathbb E[X] = \frac{r}{p} = \frac{160}{117}r \approxeq 1.3675r }\] \[\boxed{ \operatorname{Var}(X) = \frac{r(1-p)}{p^{2}} = \frac{6880}{13689}r \approxeq 0.5026r }\]袁绍（神话再临·火）
乱击：出牌阶段，你可以将任意两张相同花色的手牌当【万箭齐发】使用。
本节研究当手牌数为$n$时，可以发动乱击的次数。这里采用的方法基于二项分布。
求 $m$ 种花色的 $n$ 张手牌中，最大配对数的期望。
设 $m \ge 2$，$n \in \mathbb{N}$。令
\[X_1, X_2, \dots, X_n\]为相互独立同分布的随机变量，且
\[\mathbb{P}(X_i = j) = \frac{1}{m}, \quad j = 1,2,\dots,m\]对每一个取值$j$，定义出现次数
\[N_j := \#\{ i \in \{1,\dots,n\} : X_i = j \}\]若两个元素取值相同，则可构成一对；每个元素最多只能使用一次。因此，对取值$j$最多能形成的配对数为
\[\left\lfloor \frac{N_j}{2} \right\rfloor\]整个样本中能形成的最大配对数为
\[P_{\max} := \sum_{j=1}^{m} \left\lfloor \frac{N_j}{2} \right\rfloor\]目标是计算最大配对数的期望
\[\mathbb{E}[P_{\max}]\]虽然随机变量 $N_1,\dots,N_m$ 彼此相关，但期望具有线性性，因此
\[\mathbb{E}[P_{\max}] = \sum_{j=1}^{m} \mathbb{E}\!\left[\left\lfloor \frac{N_j}{2} \right\rfloor\right]\]由于对称性，
\[\mathbb{E}[P_{\max}] = m \, \mathbb{E}\!\left[\left\lfloor \frac{N_1}{2} \right\rfloor\right]\]注意到
\[N_1 \sim \operatorname{Bin}\!\left(n,\frac{1}{m}\right)\]对任意整数 $k$，有恒等式
\[\left\lfloor \frac{k}{2} \right\rfloor = \frac{k}{2} - \frac{1}{2}\mathbf{1}_{\{k \text{ 为奇数}\}}\]因此
\[\mathbb{E}\!\left[\left\lfloor \frac{B}{2} \right\rfloor\right] = \frac{\mathbb{E}[B]}{2} - \frac{1}{2}\mathbb{P}(B \text{ 为奇数}), \quad B \sim \operatorname{Bin}\!\left(n,\frac{1}{m}\right)\]对于 $B \sim \operatorname{Bin}(n,p)$，有
\[\mathbb{P}(B \text{ 为奇数}) = \frac{1 - (1 - 2p)^n}{2}\]在此问题中 $p = \frac{1}{m}$，因此
\[\mathbb{P}(B \text{ 为奇数}) = \frac{1 - \left(1 - \frac{2}{m}\right)^n}{2}\]由于
\[\mathbb{E}[B] = \frac{n}{m}\]可得
\[\mathbb{E}\!\left[\left\lfloor \frac{B}{2} \right\rfloor\right] = \frac{n}{2m} - \frac{1}{4} + \frac{1}{4}\left(1 - \frac{2}{m}\right)^n\]最终得到最大配对数的期望：
\[\boxed{ \mathbb{E}[P_{\max}] = \frac{n}{2} - \frac{m}{4} + \frac{m}{4}\left(1 - \frac{2}{m}\right)^n }\]当 $m=4$ 时，
\[\mathbb{E}[P_{\max}] = \frac{n}{2} - \frac{4}{4} + \frac{4}{4}\left(1 - \frac{2}{4}\right)^n = \frac{n}{2} - 1 + \frac{1}{2^{n}}\]该期望曲线为：
下面计算最大配对数的方差。
对任意整数 $u\ge 0$，有恒等式
\[\left\lfloor\frac u2\right\rfloor=\frac u2-\frac12\mathbf 1_{\{u\text{ 为奇数}\}}\]因此
\[P_{\max} =\sum_{j=1}^m\left(\frac{N_j}{2}-\frac12\mathbf 1_{\{N_j\text{ 奇}\}}\right) =\frac12\sum_{j=1}^m N_j-\frac12\sum_{j=1}^m \mathbf 1_{\{N_j\text{ 奇}\}} =\frac n2-\frac12\,O\]其中 \(O:=\sum_{j=1}^m I_j,\qquad I_j:=\mathbf 1_{\{N_j\text{ 为奇数}\}}\)
于是
\[\operatorname{Var}(P_{\max})=\frac14\,\operatorname{Var}(O)\]因此关键变成：计算 $O=\sum_{j=1}^m I_j$ 的方差。
注意到对任意整数随机变量 $Y$，有
\[\mathbf 1_{\{Y\text{ 奇}\}}=\frac{1-(-1)^Y}{2}\]所以
\[\mathbb E[I_1]=\mathbb P(N_1\text{ 奇}) =\frac{1-\mathbb E[(-1)^{N_1}]}{2}\]而 $N_1\sim\operatorname{Bin}(n,\tfrac1m)$。对于二项分布随机变量，其概率生成函数（Probability-generating function）为
\[\mathbb E[t^{N_1}]=\sum_{k=0}^{n} t^k \binom{n}{k} p^k (1-p)^{n-k} = (1-p+pt)^n,\quad p=\frac1m\]其中最后一步由二项式定理得到。
取 $t=-1$ 得到
\[\mathbb E[(-1)^{N_1}]=(1-p-p)^n=(1-2p)^n=\left(1-\frac{2}{m}\right)^n\]定义 \(a:=\left(1-\frac{2}{m}\right)^n\)，于是
\[\mathbb E[I_1]=\frac{1-a}{2}\]由于 $I_1$ 是 $0/1$ 变量，有 $I_1^2 = I_1$，因此$\mathbb{E}[I_1^2] = \mathbb{E}[I_1]$
所以
\[\begin{aligned} \operatorname{Var}(I_1) &= \mathbb{E}[I_1^2] - (\mathbb{E}[I_1])^2 \\ &= \mathbb E[I_1](1-\mathbb E[I_1]) \\ &= \frac{1-a}{2}\cdot\frac{1+a}{2} \\ &= \frac{1-a^2}{4} \end{aligned}\]计算协方差 $\operatorname{Cov}(I_1,I_2)$，需要 $\mathbb E[I_1I_2]$：
\[I_j=\frac{1-(-1)^{N_j}}{2} \quad\Longrightarrow\quad I_1I_2=\frac{1-(-1)^{N_1}-(-1)^{N_2}+(-1)^{N_1+N_2}}{4}\]取期望：
\[\mathbb E[I_1I_2] =\frac{1-\mathbb E[(-1)^{N_1}]-\mathbb E[(-1)^{N_2}]+\mathbb E[(-1)^{N_1+N_2}]}{4}\]由对称性 $\mathbb E[(-1)^{N_1}]=\mathbb E[(-1)^{N_2}]=a$，所以
\[\mathbb E[I_1I_2]=\frac{1-2a+\mathbb E[(-1)^{N_1+N_2}]}{4}\]现在计算 $\mathbb E[(-1)^{N_1+N_2}]$。令 \(S:=N_1+N_2.\)
由于每张牌的花色独立地以概率 $\frac{1}{m}$ 取每个值，故花色落在 \(\{1,2\}\) 的概率为 $\frac{2}{m}$。因此
\[S = N_1 + N_2 \sim\operatorname{Bin}\!\left(n,\frac{2}{m}\right)\]利用前面的概率生成函数结果（取 $p = \frac{2}{m}$，$t = -1$）：
\[\mathbb E[(-1)^S] =\left(1-2 \cdot \frac{2}{m}\right)^n =\left(1-\frac{4}{m}\right)^n\]定义 \(b:=\left(1-\frac{4}{m}\right)^n\) ，代回得到
\[\mathbb E[I_1I_2]=\frac{1-2a+b}{4}\]于是
\[\operatorname{Cov}(I_1,I_2) =\mathbb E[I_1I_2]-\mathbb E[I_1]^2 =\frac{1-2a+b}{4}-\left(\frac{1-a}{2}\right)^2 =\frac{b-a^2}{4}\]因为 $O=\sum_{j=1}^m I_j$，且各 $I_j$ 交换对称，所以
\[\begin{aligned} \operatorname{Var}(O) &= \sum_{j=1}^m \operatorname{Var}(I_j) + \sum_{j \neq k} \operatorname{Cov}(I_j, I_k) \\ &= m\,\operatorname{Var}(I_1)+m(m-1)\,\operatorname{Cov}(I_1,I_2) \end{aligned}\]代入上面的结果：
\[\operatorname{Var}(O) = m\cdot\frac{1-a^2}{4}+m(m-1)\cdot\frac{b-a^2}{4} =\frac{m}{4}(1-a^2)+\frac{m(m-1)}{4}(b-a^2)\]因此
\[\boxed{ \operatorname{Var}(P_{\max}) =\frac14\,\operatorname{Var}(O) = \frac{m}{16}(1-a^2)+\frac{m(m-1)}{16}(b-a^2) }\]其中
\[a=\left(1-\frac{2}{m}\right)^n,\qquad b=\left(1-\frac{4}{m}\right)^n\]当 $m=4$ 时
\[a=\left(\frac12\right)^n,\qquad b=(1-1)^n=0\quad(n\ge 1)\]代入得到
\[\operatorname{Var}(P_{\max}) =\frac{4}{16}\left(1-2^{-2n}\right)+\frac{12}{16}\left(0-2^{-2n}\right) =\frac{1}{4}-\frac{1}{4^n}\]对于袁绍的技能乱击，我们还会研究另一个课题。见“生日问题”一节）。
二次函数
王朗（SP）
鼓舌：出牌阶段限一次，你可以用一张手牌与至多三名角色同时拼点，然后依次结算拼点结果，没赢的角色选择一项：1. 弃置一张牌；2. 令你摸一张牌。若拼点没赢的角色是你，你需先获得一个“饶舌”标记（你有7个“饶舌”标记时，你死亡）。
激词：当你发动“鼓舌”拼点的牌亮出后，若点数小于X，你可令点数+X；若点数等于X，你可以令你本回合发动“鼓舌”的次数上限+1。（X为你“饶舌”标记的数量）
下面考察发动激词的情况下拼点赢的概率。定义王朗玩家的拼点牌点数为$T$，“饶舌”标记数量为$x$，则拼点赢的可能情况数为：
\[C_{T,x} = \min\{13,\; T + x\mathbf{1}_{\{T<x\}}\} - 1\]拼点赢的概率：
\[\begin{aligned} p(x) &= \mathbb E_T\!\left[\mathbb P\!\left(U < T + \mathbf{1}_{\{T<x\}}x \mid T\right)\right] \\ &= \frac{1}{13}\sum_{T=1}^{13}\frac{C_{T,x}}{13} \\ &= \frac{1}{169}\left(\sum_{i=1}^{x-1}(i+x-1)+\sum_{i=x}^{13}(i-1)\right) \\ &= \frac{x^2-x+78}{169} \end{aligned}\]其中$U$表示对手的拼点牌点数，$0 \le x \le 6$。
由于对$T<x$的区间求和产生了等差数列求和项，故$p(x)$关于$x$为二次函数。
其概率散点图如下：
可见，随着“饶舌”标记的增多，拼点赢的概率连续增大。
接下来研究“饶舌”标记数量为$x$时，发动鼓舌而获得的“饶舌”标记数量$K$。假设发动鼓舌时与$n$名角色同时拼点，其中$0 < n \le 3$。很明显，其概率符合二项分布$K \sim \operatorname{Bin}(n,1-p(x))$，数学期望和方差分别为：
\[\boxed{\mathbb E[K;x] = n(1-p(x)) = \frac{n}{169}(91-x^2+x)}\] \[\boxed{ \begin{aligned} \operatorname{Var}(K;x) &= np(x)(1-p(x)) \\ &= \frac{n}{28561}(91-x^2+x)(x^2-x+78) \\ &= \frac{n}{28561}(-x^4+2x^3+12x^2-13x+7098) \end{aligned} }\]其期望散点图如下：
以$n = 3$为例：
最后，计算发动鼓舌后角色死亡的概率：
\[\mathbb{P}(K \ge 7-x) = 1-\mathbb{P}(K \le 6-x) = 1-F(6-x;n;1-p(x)) = 1-\sum_{i=0}^{6-x}\binom{n}{i}(1-p(x))^{i}p(x)^{n-i}\]其中$F$为累积分布函数（Cumulative distribution function）。
\[\begin{aligned} \mathbb{P}(K \ge 7-x) = \begin{cases} 0 & \text{if $x \le 6-n$} \\ 1-\frac{1}{169^{n}}\sum_{i=0}^{6-x}\binom{n}{i}(91-x^2+x)^{i}(x^2-x+78)^{n-i} & \text{otherwise} \end{cases} \end{aligned}\]其概率散点图如下：
以$n = 3$为例：
生日问题
孙资&刘放（原创之魂2016）
（图片来源：https://www.bilibili.com/read/cv16369119）
瑰藻：若你于弃牌阶段弃置过至少两张牌且这些牌花色均不相同，则你可以回复1点体力或摸一张牌。
下面我们来计算瑰藻可以回复体力或摸牌的概率。
定义弃牌阶段弃置的牌数为$n \enspace (n > 1)$，一副牌中不同的花色数为$m$。在理想模型下，可将弃置牌的花色近似视为独立均匀分布于$m$类。此问题为典型的生日问题（Birthday problem）：
\[{\bar {p}}(n) = \frac{m!}{(m-n)!m^{n}}, \quad 1 < n \le m\]根据鸽巢原理（Pigeonhole principle），当$n > m$时，${\bar {p}}(n) = 0$。
在此模型下，$m=4$。
其概率散点图为：
周泰（神话再临·风）
（图片来源：https://www.bilibili.com/read/cv15635764）
不屈：锁定技，当你处于濒死状态时，你将牌堆顶的一张牌置于你的武将牌上，称为“创”，若此牌的点数与已有的“创”点数均不同，则你将体力回复至1点。若出现相同点数则将此牌置入弃牌堆，若你的武将牌上有“创”，则你的手牌上限与“创”的数量相等。
我们考虑不屈失败时的“创”数。定义$X$为出现首次点数重复所需翻牌次数（包含触发重复的那一次），一副牌中不同的点数为$M$。此问题可用生日问题的一种形式来描述。
PMF：
\[p_{X}(k) = \frac{m!(k-1)}{m^{k}(m-k+1)!}, \quad k \in \mathbb{Z} \enspace \text{and} \enspace 1 \le k \le m + 1\]在此模型下，$m=13$。
其概率散点图为：
数学期望为：
\[\boxed{ \mathbb E[X] = 1 + \sum_{k=1}^{m}{\frac{m!}{(m-k)!m^{k}}} \approxeq 5.2124 }\]我们也可以将右式进行拉马努金渐近展开，从而得到简化的近似表达式：
\[\mathbb E[X] \sim \sqrt{\frac{\pi m}{2}}+\frac {2}{3} \approxeq 5.1856\]由此我们可知，周泰的“创”数平均约为5.2张。
接下来我们考察周泰的Online战功“金枪不倒”：使用周泰在一局游戏中拥有过至少9张“不屈”牌并且未死。
\[p_{X}(k \ge 9) = 1 - \sum_{k = 2}^8{p_X(k)} \approxeq 0.06361\]【点数非均匀修正】
若考虑【2】与【Q】的概率较于其他点数略高，我们可以对PMF做细微修正，得到更加准确的数学期望：
\[\mathbb E[X] \approxeq 5.2065\]可见，此结果的差异完全可以忽略。
神吕蒙（神话再临·风）
涉猎：摸牌阶段，你可选择采取以下行动来取代摸牌：从牌堆顶亮出五张牌，拿走不同花色的牌各一张，弃掉其余的。
下面我们研究涉猎获得的牌数$X$。该问题是生日问题的一种形式。根据第二类斯特林数（Stirling numbers of the second kind）的定义：
PMF：
\[p_{X}(k) = \frac{m^{\underline k}}{m^{n}}\begin{Bmatrix}n\\k\end{Bmatrix}\]其中$n$为亮出的牌数，$m$为一副牌中不同的花色数。在此模型下，$n=5, \enspace m=4$。
其概率散点图为：
数学期望和方差分别为：
\[\boxed{ \mathbb E[X] = m(1-(1-\frac{1}{m})^n) = \frac{781}{256} \approxeq 3.0508 }\] \[\boxed{ \operatorname{Var}(X) = m(m-1)(1-\frac{2}{m})^{n}+m(1-\frac{1}{m})^{n}-m^2(1-\frac{1}{m})^{2n} = \frac{27735}{65536} \approxeq 0.4232 }\]袁绍（神话再临·火）
对于袁绍的技能乱击，除了前面提到的最大配对数，我们还可以研究发动乱击的概率。
定义手牌数为$X$，一副牌中不同的花色数为$m$。可以发动乱击的概率可用生日问题的另一种形式来表示：
\[p_{X}(n) = 1-\prod_{k=1}^{n-1}\left(1-{\frac{\min(k, m)}{m}}\right)\]在此模型下，$m=4$。
其概率散点图为：
和约束下可选子集最大规模的期望
曹冲（新一将成名2013）
（图片来源：https://www.bilibili.com/read/cv15576366）
称象：当你受到伤害后，你可以亮出牌堆顶的四张牌，然后你获得其中的任意张点数之和不大于13的牌。
下面我们考虑理想模型下称象获得的牌数。
【Ehrhart】
用 \(\mathbf{x} = (x_1,x_2,x_3,x_4)^{\mathsf T},\qquad x_i \in \{\,1,2,\dots\,\}\) 表示亮出的四张牌的点数。定义随机变量$K$为称象最终获得的牌数，即满足点数和不超过 13 的可选牌张数。
【情形一】
当$K = 4$时，四张牌全部可取，需满足
\[x_1+x_2+x_3+x_4 \le 13\]运用Stars and bars法，可知其整数解组数为
\[\binom{13}{4} = 715\]【情形二】
当$K = 1$时，当且仅当恰有一张牌可取，对应的不等式系统可写为
\[A\mathbf x\le \mathbf b\]其中
\[A=\begin{pmatrix} -1 & -1 & 0 & 0 \\ -1 & 0 & -1 & 0 \\ -1 & 0 & 0 & -1 \\ 0 & -1 & -1 & 0 \\ 0 & -1 & 0 & -1 \\ 0 & 0 & -1 & -1 \\ 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}, \qquad \mathbf b= \begin{pmatrix} -14\\-14\\-14\\-14\\-14\\-14\\ 13\\13\\13\\13\\ -1\\-1\\-1\\-1 \end{pmatrix}\]该不等式组定义了一个四维凸多胞形（Convex polytope）。问题等价于计算该多胞形中晶格点的数量。以上矩阵不等式为该四维凸多胞形的半空间表示（Half-space representation）。利用开源软件polymake计算可知，该凸多胞形为整多胞形（Integral polytope），其Ehrhart多项式（Ehrhart polynomial）：
1
2
3
4
5
6
7
polytope > $p = new Polytope(INEQUALITIES=>[[-14,1,1,0,0],[-14,1,0,1,0],[-14,1,0,0,1],[-14,0,1,1,0],[-14,0,1,0,1],[-14,0,0,1,1],[-1,1,0,0,0],[-1,0,1,0,0],[-1,0,0,1,0],[-1,0,0,0,1],[13,-1,0,0,0],[13,0,-1,0,0],[13,0,0,-1,0],[13,0,0,0,-1]]);
polytope > print $p->LATTICE;
true
polytope > print $p->EHRHART_POLYNOMIAL;
2592*x^4 + 1296*x^3 + 252*x^2 + 24*x + 1
因此，该凸多胞形内晶格点的个数为
\[\#(K=1)=2592+1296+252+24+1=4165\]【情形三】
当$K = 3$时，根据不可取牌在四个位置中的不同分布情况，可将可行域划分为若干互不相交的情形。每一种情形均可表示为一个线性不等式系统
\[A_i \mathbf x \le \mathbf b_i\]从而对应一个四维凸多胞形。
\[A_1 = \begin{pmatrix} -1 & -1 & -1 & -1 \\ 1 & 1 & 1 & 0 \\ 1 & 1 & 0 & 1 \\ 1 & 0 & 1 & 1 \\ 0 & 1 & 1 & 1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}, \qquad \mathbf b_1 = \begin{pmatrix} -14 \\ 13 \\ 13 \\ 13 \\ 13 \\ -1 \\ -1 \\ -1 \\ -1 \end{pmatrix}\] \[A_2 = \begin{pmatrix} 1 & 1 & 1 & 0 \\ 1 & 1 & 0 & 1 \\ 1 & 0 & 1 & 1 \\ 0 & -1 & -1 & -1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}, \qquad \mathbf b_2 = \begin{pmatrix} 13 \\ 13 \\ 13 \\ -14 \\ -1 \\ -1 \\ -1 \\ -1 \end{pmatrix}\] \[A_3 = \begin{pmatrix} 1 & 1 & 1 & 0 \\ 1 & 1 & 0 & 1 \\ -1 & 0 & -1 & -1 \\ 0 & -1 & -1 & -1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}, \qquad \mathbf b_3 = \begin{pmatrix} 13 \\ 13 \\ -14 \\ -14 \\ -1 \\ -1 \\ -1 \\ -1 \end{pmatrix}\] \[A_4 = \begin{pmatrix} -1 & -1 & -1 & -1 \\ 1 & 1 & 1 & 0 \\ -1 & -1 & 0 & -1 \\ -1 & 0 & -1 & -1 \\ 0 & -1 & -1 & -1 \\ 0 & 0 & 0 & 1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}, \qquad \mathbf b_4 = \begin{pmatrix} -14 \\ 13 \\ -14 \\ -14 \\ -14 \\ 13 \\ -1 \\ -1 \\ -1 \\ -1 \end{pmatrix}\]上述四个不等式系统分别给出了四个四维凸多胞形。
由牌位的对称性可知，其对应的组合系数分别为：$\binom{4}{4} = 1$、$\binom{4}{3} = 4$、$\binom{4}{2} = 6$和$\binom{4}{1} = 4$。
再次利用polymake的计算可知，以上四种情况对应的四个凸多胞形均为有理多胞形（Rational polytope）而非整多胞形。因此，只能通过求解其Ehrhart准多项式（Ehrhart quasi-polynomial）来求解其所包含的晶格点：
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
polytope > $p1 = new Polytope(INEQUALITIES=>[[-14,1,1,1,1],[13,-1,-1,-1,0],[13,-1,-1,0,-1],[13,-1,0,-1,-1],[13,0,-1,-1,-1],[-1,1,0,0,0],[-1,0,1,0,0],[-1,0,0,1,0],[-1,0,0,0,1]]);
polytope > $p2 = new Polytope(INEQUALITIES=>[[13,-1,-1,-1,0],[13,-1,-1,0,-1],[13,-1,0,-1,-1],[-14,0,1,1,1],[-1,1,0,0,0],[-1,0,1,0,0],[-1,0,0,1,0],[-1,0,0,0,1]]);
polytope > $p3 = new Polytope(INEQUALITIES=>[[13,-1,-1,-1,0],[13,-1,-1,0,-1],[-14,1,0,1,1],[-14,0,1,1,1],[-1,1,0,0,0],[-1,0,1,0,0],[-1,0,0,1,0],[-1,0,0,0,1]]);
polytope > $p4 = new Polytope(INEQUALITIES=>[[13,-1,-1,-1,0],[-14,1,1,0,1],[-14,1,0,1,1],[-14,0,1,1,1],[-1,1,0,0,0],[-1,0,1,0,0],[-1,0,0,1,0],[-1,0,0,0,1],[13,0,0,0,-1]]);
polytope > print $p1->LATTICE;
false
polytope > print $p2->LATTICE;
false
polytope > print $p3->LATTICE;
false
polytope > print $p4->LATTICE;
false
polytope > print $p1->EHRHART_QUASI_POLYNOMIAL;
1250/9*x^4 + 1750/9*x^3 + 1675/18*x^2 + 35/2*x + 1
1250/9*x^4 + 1750/9*x^3 + 1675/18*x^2 + 35/2*x + 10/9
1250/9*x^4 + 1750/9*x^3 + 1675/18*x^2 + 35/2*x + 1
polytope > print $p2->EHRHART_QUASI_POLYNOMIAL;
256/9*x^4 + 512/9*x^3 + 352/9*x^2 + 32/3*x + 1
256/9*x^4 + 512/9*x^3 + 352/9*x^2 + 32/3*x + 8/9
256/9*x^4 + 512/9*x^3 + 352/9*x^2 + 32/3*x + 1
polytope > print $p3->EHRHART_QUASI_POLYNOMIAL;
729/8*x^4 + 243/2*x^3 + 117/2*x^2 + 12*x + 1
729/8*x^4 + 243/2*x^3 + 117/2*x^2 + 12*x + 7/8
polytope > print $p4->EHRHART_QUASI_POLYNOMIAL;
5875/9*x^4 + 4400/9*x^3 + 1240/9*x^2 + 53/3*x + 1
5875/9*x^4 + 4400/9*x^3 + 1240/9*x^2 + 53/3*x + 8/9
5875/9*x^4 + 4400/9*x^3 + 1240/9*x^2 + 53/3*x + 1
polytope > print $p1->N_LATTICE_POINTS;
445
polytope > print $p2->N_LATTICE_POINTS;
136
polytope > print $p3->N_LATTICE_POINTS;
284
polytope > print $p4->N_LATTICE_POINTS;
1298
在乘以相应的系数后，此类整数解组数为
\[\#(K=3) = 445+136*4+284*6+1298*4=7885\]【情形四】
当$K = 2$时，由于所有点数组合的总数为 (13^4)，可通过补集直接得到
\[\#(K=2)=13^4-\#(K=1)-\#(K=3)-\#(K=4)=15796\]由上述计数结果可得随机变量$K$的概率分布。其概率散点图为：
数学期望和方差分别为：
\[\boxed{\mathbb E[K] = \frac{62272}{28561} \approxeq 2.1803}\] \[\boxed{\operatorname{Var}(K) = \frac{399322010}{815730721} \approxeq 0.4895}\]该问题亦可通过组合生成函数（Generating function）与次序统计量进行纯组合论推导（参见附录）。但在本问题规模下，该方法并不能显著降低计算复杂度；相反，将计数问题转化为凸多胞形中的晶格点计数，并借助 Ehrhart 理论与 polymake 实现自动化计算，更有利于保证结果的完整性与可复现性。
随机子集和的可达性概率
糜竺（SP）
资援：出牌阶段限一次，你可以将任意张点数之和为13的手牌交给一名其他角色，然后该角色回复1点体力。
下面研究在手牌数为 $X$ 的情况下，可以发动资援的概率。这是典型的NP完备（NP-Complete）的子集和问题（Subset sum problem）。程序如下：
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
from collections import defaultdict
from math import gcd
def subset_sum_probability(n, target=13, max_value=13):
"""
Computes the exact probability that among n i.i.d. draws from {1,2,...,max_value},
there exists a nonempty subset whose sum equals `target`.
The result is returned as a reduced fraction (numerator, denominator).
Core idea:
- Track which subset sums (0..target) are achievable so far using a bitmask.
- Use dynamic programming over masks.
"""
# We only care about subset sums from 0 to `target`
# Bit i in the mask means: "sum i is achievable"
WIDTH = target + 1 # number of bits we care about
FULL_MASK = (1 << WIDTH) - 1 # mask with bits 0..target set to 1
TARGET_BIT = 1 << target # bit representing sum == target
# ------------------------------------------------------------------
# Precompute state transitions
#
# next_mask[mask][x] = new mask after adding value x to the multiset
#
# Transition rule:
# new_sums = old_sums ∪ (old_sums + x)
#
# In bit operations:
# new_mask = mask | (mask << x)
# and then we truncate to keep only bits <= target
# ------------------------------------------------------------------
next_mask = [[0] * (max_value + 1) for _ in range(FULL_MASK + 1)]
for mask in range(FULL_MASK + 1):
for x in range(1, max_value + 1):
shifted = (mask << x) & FULL_MASK
next_mask[mask][x] = mask | shifted
# ------------------------------------------------------------------
# DP array:
# dp[mask] = number of ordered sequences (x1,...,xk)
# that produce exactly this subset-sum mask
#
# Initially, before seeing any numbers:
# - Only sum 0 is achievable
# - mask = 1 << 0 = 1
# ------------------------------------------------------------------
dp = [0] * (FULL_MASK + 1)
dp[1] = 1
# ------------------------------------------------------------------
# Main DP loop: add one draw at a time
# Each step corresponds to drawing one more number Xi ∈ {1..max_value}
# ------------------------------------------------------------------
for _ in range(n):
new_dp = [0] * (FULL_MASK + 1)
# For each current state (mask) and its count
for mask, count in enumerate(dp):
if count == 0:
continue
# Try all possible next values x
row = next_mask[mask]
for x in range(1, max_value + 1):
new_mask = row[x]
new_dp[new_mask] += count
dp = new_dp
# ------------------------------------------------------------------
# Count successful sequences:
# those whose final mask has the TARGET_BIT set,
# meaning "there exists a subset summing to target"
# ------------------------------------------------------------------
numerator = sum(
count for mask, count in enumerate(dp)
if mask & TARGET_BIT
)
# Total number of possible sequences
denominator = max_value ** n
# Reduce fraction
g = gcd(numerator, denominator)
return numerator // g, denominator // g
该程序的时间复杂度为 $\mathcal{O}(n \cdot 2^{T+1} \cdot m)$，空间复杂度为 $\mathcal{O}(2^{T+1} \cdot m)$。其中 $T$ 为目标和（$13$），$m$为卡牌最大值（$13$）。由于计算量大，我们只计算 $0 < X \le 10$ 的情况：
1
2
3
4
5
6
7
8
9
10
1
37
931
18485
303026
4391122
60044006
799241127
10505121928
137264030196
其概率散点图为：
【封闭形式】
下面来说明该概率必然存在封闭形式（Closed-form Expression）.
我们研究的是概率
\[p_n = \mathbb{P}\Big( \exists\, I\subseteq\{1,\dots,n\},\ I\neq\varnothing: \sum_{i\in I} X_i = 13 \Big)\]其中, $X_i$ 独立同分布，取值于 ${1,\dots,13}$
核心问题是：是否存在某个子集和等于 13。这个“是否存在”的判定，并不依赖于所有历史细节，而只依赖于一件事：当前已经“可达”的子集和有哪些（在 0 到 13 之间）
定义状态：
\[R_n := \{\text{前 } n \text{ 个数能组成的所有子集和}\}\cap\{0,1,\dots,13\}\]每个状态可以用一个 14 位二进制掩码（bitmask）表示。第 $k$ 位为 1 表示“子集和 $k$ 是可达的”。那么，总状态数最多为：
\[2^{14} = 16384\]即，状态空间是一个有限集合。
当加入第 $n+1$ 个随机变量 $X_{n+1}=x$ 时，状态更新规则是：
\[R_{n+1} = R_n \;\cup\; \{ r+x : r\in R_n,\ r+x\le 13 \}\]该规则只依赖于 $s$ ，与 $n$ 无关。其随机性只来自输入 $x$，且 $x$ 的分布在每一步都是相同的（均匀分布）。
因此，这是一个有限状态、时间齐次的马尔可夫链（Time-homogeneous Markov chain with a finite state space）。进而我们知道，该概率值必然存在封闭形式。
下面举两个封闭形式的例子：
一、有限维线性递推
由于状态演化可以写成
\[\mathbf{v}_{n+1} = \mathbf{v}_n P\]其中，$\mathbf{v}_n$ 为状态分布向量（有限维）， $P$ 为固定的转移矩阵。
$\mathbf{v}_n$ 的任意线性函数（例如 $p_n$ ） 都满足一个有限阶、常系数的线性递推关系
即存在 $L<\infty$ 和常数 $c_1,\dots,c_L$ ，使得
\[p_n = c_1 p_{n-1} + \cdots + c_L p_{n-L} \quad (n\ \text{充分大})\]二、生成函数是有理函数
定义生成函数： \(G(z) = \sum_{n\ge 0} p_n z^n\)
由于系统是有限维线性的，可以证明：
\(G(z) = \frac{A(z)}{B(z)}\) 其中 $A,B$ 是多项式。
【证明】
下面证明：可以发动资援的概率随着手牌数单调递增，且其极限为$1$。
命题1： 令
\[E_n := \Big\{\exists\, I\subseteq\{1,2,\dots,n\},\ I\neq\varnothing:\ \sum_{i\in I}X_i = 13\Big\}, \qquad p_n := \mathbb{P}(E_n)\]其中 $X_1,X_2,\dots$ 为独立同分布随机变量，且
\[X_i \sim \mathrm{Unif}\{1,2,\dots,13\}\]则序列 $(p_n)_{n\ge 1}$ 严格单调递增，即
\[p_{n+1} > p_n \quad \text{对所有 } n\ge 1 \text{成立}\]证明：
首先注意到对任意 $n$，都有集合包含关系
\[E_n \subseteq E_{n+1}\]因此由概率单调性原理或测度单调性（Measure Monotonicity），
\[p_n \le p_{n+1}\]下面证明不等号是严格的。
考虑如下事件：
\[B_n := \Big\{ \sum_{i\in I} X_i \neq 13 \text{ 对所有 } I\subseteq\{1,\dots,n\},\ I\neq\varnothing, \ \text{但 } X_{n+1}=13 \Big\}\]在事件 $B_n$ 中：
- 前 $n$ 个随机变量无法组成和为 $13$ 的子集，因此 $E_n$ 不发生；
- 但由于 $X_{n+1}=13$，取单元素集合 ${n+1}$ 即可得到和为 $13$ ，因此 $E_{n+1}$ 发生。
故
\[B_n \subseteq E_{n+1}\setminus E_n\]由于
\[\mathbb{P}(B_n) = \mathbb{P}\big(X_{n+1}=13\big)\, \mathbb{P}\big(E_n^c\big) = \frac{1}{13}\,\mathbb{P}(E_n^c) > 0\]（注意 $\mathbb{P}(E_n^c)>0$，例如当 $X_1=\cdots=X_n=2$ 时，任意非空子集和都是偶数，不可能等于 $13$，所以 $E_n^c$ 必然发生），
因此
\[\mathbb{P}(E_{n+1}) - \mathbb{P}(E_n) = \mathbb{P}(E_{n+1}\setminus E_n) \ge \mathbb{P}(B_n) > 0\]从而
\[p_{n+1} > p_n\] \[\square\]命题2： \(\lim_{n\to\infty} p_n = 1\)
证明：
定义事件
\[A_n := \{\exists\, i\in\{1,\dots,n\}: X_i = 13\}\]若 $A_n$ 发生，则取单元素集合 $I={i}$ 就有 $\sum_{j\in I}X_j = X_i = 13$，因此
\[A_n \subseteq E_n\]于是
\[p_n = \mathbb{P}(E_n) \ge \mathbb{P}(A_n)\]接下来计算 $\mathbb{P}(A_n)$。由于各 $X_i$ 独立且 $\mathbb{P}(X_i\neq 13)=\frac{12}{13}$，有
\[\mathbb{P}(A_n) = 1-\mathbb{P}(\forall i,\ X_i\neq 13) = 1-\left(\frac{12}{13}\right)^n\]另一方面显然 $p_n\le 1$（因为是概率），因此得到夹逼不等式
\[1-\left(\frac{12}{13}\right)^n \le p_n \le 1\]当 $n\to\infty$ 时， $\left(\frac{12}{13}\right)^n \to 0$，从而左端趋于 $1$。由夹逼定理（Squeeze Theorem）可得
\[\lim_{n\to\infty} p_n = 1\] \[\square\]从另一个角度看，可达子集和集合的演化构成一个有限状态、时间齐次的吸收型马尔可夫链（Absorbing Markov Chain）。所有包含目标值 $13$ 的状态都是吸收态，一旦进入便无法离开。由于从初始状态出发进入吸收态的概率为正，且不存在非吸收的封闭类，该链以概率 $1$ 最终被吸收。因此，目标事件的概率 $p_n$ 随 $n$ 增大必然趋于 $1$。
技能比较
最后，我们在理想模型下比较一些类似的技能。
额外摸牌
| 英姿 | 闭月 | 洛神 | 涉猎 | 再起 | |
| 期望 | 1 | 1 | 1 | 1.0508 | 0.25X |
| 方差 | 0 | 0 | 2 | 0.4232 | 0.1875X |
注1：再起中的$X$表示已损失体力值。 注2：方差越大表示收益波动越大，越不稳定。
- 英姿和闭月带来稳定的额外摸牌收益
- 洛神、涉猎和再起的额外摸牌不稳定
- 洛神、涉猎和再起的额外摸牌受到花色条件的约束
- 涉猎和再起需要亮出所得的牌
- 再起的摸牌数一般小于1，但有一定概率可以回复体力
受到伤害后摸牌
| 遗计 | 称象 | |
| 期望 | 2 | 2.1803 |
| 方差 | 0 | 0.4895 |
- 遗计的摸牌数稳定
- 称象的期望略高，但摸牌数不稳定
- 称象的摸牌有点数的限制
总结
以上几个同类技能各有优劣，可见设计是比较均衡的。
附录
称象获得牌数的组合生成函数解法
另一种解法从组合生成函数出发，将“无序多重集合计数”与“有序四元组计数”通过一个明确的升标号算子连接起来。该方法完全白盒，且所有计算均归结为有限次系数抽取与有限求和。
设四张牌点数为
\[x_1,x_2,x_3,x_4\in\{1,2,\dots,13\}\]并记其排序后的次序统计量为
\[A_1\le A_2\le A_3\le A_4\]定义在总和不超过 $13$ 的约束下最多可取的张数
\[K:=\max\{k\in\{1,2,3,4\}:A_1+\cdots+A_k\le 13\}\]显然
\[\mathbb P(K\ge k)=\mathbb P(A_1+\cdots+A_k\le 13)\]并且
\[\mathbb P(K=k)=\mathbb P(K\ge k)-\mathbb P(K\ge k+1)\]因此，只需分别计数事件 $K\ge 2$、$K\ge 3$、$K\ge 4$ 对应的有序四元组数量。
多重集合的单项式表示
引入形式变量 $y_1,\dots,y_{13}$，用于标记各点数的出现情况。
对任意大小为 $4$ 的多重集合
\[M=\{1^{m_1}2^{m_2}\cdots 13^{m_{13}}\},\qquad \sum_{j=1}^{13}m_j=4\]其中 $m_j$ 表示点数 $j$ 在多重集合中的重数（出现次数）。定义其对应单项式
\[y^M:=\prod_{j=1}^{13}y_j^{m_j}\]例：多重集合 ${1,1,2,5}$ 对应 $m_1=2, m_2=1, m_5=1$（其余 $m_j=0$），其单项式为 $y_1^2 y_2 y_5$。
升标号算子
升标号（Labeling）是指将无序的多重集合分配到有序的位置上。
定义升标号线性算子 $\mathrm{Lab}_4$，其在基单项式上的作用为
\[\mathrm{Lab}_4\!\left(\prod_{j=1}^{13}y_j^{m_j}\right) :=\frac{4!}{\prod_{j=1}^{13}m_j!}\prod_{j=1}^{13}y_j^{m_j}, \qquad \sum m_j=4\]并线性延拓。
组合意义：系数 $\dfrac{4!}{\prod m_j!}$ 正是多项式系数（multinomial coefficient）
\[\binom{4}{m_1, m_2, \dots, m_{13}} = \frac{4!}{\prod_{j=1}^{13} m_j!}\]它计算的是将该多重集合的 $4$ 个元素分配到 $4$ 个有区分的位置上的排列数。
例：多重集合 ${1,1,2}$ 的升标号系数为 $\dfrac{3!}{2!\cdot 1!}=3$，对应的 $3$ 个有序排列为 $(1,1,2), (1,2,1), (2,1,1)$。
从多重集合到有序四元组的计数
对于任一仅依赖多重集合（即对置换对称）的性质 $\mathcal P$，定义多重集合枚举多项式
\[U_{\mathcal P}(y_1,\dots,y_{13}) :=\sum_{\substack{M:\ \sum m_j=4\\ \mathcal P(M)\ \text{成立}}} y^M\]则满足性质 $\mathcal P$ 的有序四元组数量为
\[\#\{(x_1,x_2,x_3,x_4):\mathcal P\} = \Big(\mathrm{Lab}_4(U_{\mathcal P})\Big)\Big|_{y_1=\cdots=y_{13}=1}\]构造 $K\ge k$ 的多重集合枚举多项式
按第 $k$ 小值分层
事件 $K\ge k$ 等价于 $A_1+A_2+\cdots+A_k\le 13$，这仅依赖多重集合结构。
按第 $k$ 小值 $A_k=c$ 进行分层。对于固定的 $c\in{1,2,\dots,13}$：
- 前 $k$ 张牌满足 $A_1\le A_2\le \cdots \le A_k = c$，即都取自 ${1,\dots,c}$ 且最大值恰为 $c$
- 第 $k$ 张牌之后的牌（若有）满足 $\ge c$
生成函数构造
为了用生成函数同时记录张数与点数和，引入标记变量：
- $u$：记录取了多少张牌
- $z$：记录点数之和
对每个点数 $j$，若它在多重集合中出现 $t$ 次，则贡献单项式 $(u\,y_j\,z^j)^t$。
点数 $j$ 可出现任意多次（包括 $0$ 次）的生成函数为几何级数：
\[\sum_{t=0}^{\infty}(u\,y_j\,z^j)^t = \frac{1}{1-u\,y_j\,z^j}\]点数 $j$ 至少出现一次的生成函数为：
\[\sum_{t=1}^{\infty}(u\,y_j\,z^j)^t = \frac{u\,y_j\,z^j}{1-u\,y_j\,z^j}\]$\Theta_k(c)$ 的定义
定义辅助生成函数
\[B_{c-1}(u,z):=\prod_{j=1}^{c-1}\frac{1}{1-u\,y_j\,z^j}\]它枚举所有点数严格小于 $c$ 的牌的任意组合。
定义截断算子 $[z^{\le 13}]$，表示只保留 $z$ 的幂次不超过 $13$ 的项：
\[[z^{\le 13}]\Phi(z):=\sum_{s=0}^{13}[z^s]\Phi(z)\]则按 $A_k=c$ 分层的前 $k$ 张牌的多重集合枚举多项式为
\[\Theta_k(c) := [z^{\le 13}]\,[u^k]\; \left(\frac{u\,y_c\,z^c}{1-u\,y_c\,z^c}\right) B_{c-1}(u,z)\]解释：
- $\dfrac{u\,y_c\,z^c}{1-u\,y_c\,z^c}$：点数 $c$ 至少出现一次（保证 $A_k=c$）
- $B_{c-1}(u,z)$：点数 $1,\dots,c-1$ 各出现任意多次
- $[u^k]$：恰好取 $k$ 张牌
- $[z^{\le 13}]$：前 $k$ 张牌的点数和不超过 $13$
$K\ge k$ 的完整公式
剩余的 $4-k$ 张牌（即 $A_{k+1}, \dots, A_4$）必须满足 $\ge A_k = c$，即取自 ${c, c+1, \dots, 13}$。
因此，完整的多重集合枚举多项式为
\[U_{\ge k}(y) = \sum_{c=1}^{13}\Theta_k(c) \cdot \left(\sum_{d=c}^{13}y_d\right)^{4-k}\]对应的有序四元组计数为
\[\#(K\ge k) = \Big(\mathrm{Lab}_4(U_{\ge k})\Big)\Big|_{y_1=\cdots=y_{13}=1}\]代码及计算结果
代码：
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
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
"""
Compute #(K >= k) using generating functions directly.
Key formula:
Theta_k(c) = [z^{<=13}][u^k] * (u*y_c*z^c / (1-u*y_c*z^c)) * prod_{j=1}^{c-1} 1/(1-u*y_j*z^j)
The y_j variables track multiplicities, enabling:
1. Decoding of concrete multisets from coefficients
2. Direct computation of labeling coefficient Lab_n = n! / prod(m_j!)
"""
from collections import defaultdict
from math import factorial
from itertools import combinations_with_replacement
def multiply_poly(p1, p2, c, max_u, max_z):
"""
Multiply two polynomials p1 * p2.
Polynomial representation: {(u_power, z_power, y_tuple): coefficient}
where y_tuple = (m_1, ..., m_c) records multiplicity of each card value.
"""
result = defaultdict(int)
for (u1, z1, y1), c1 in p1.items():
for (u2, z2, y2), c2 in p2.items():
u_new = u1 + u2
z_new = z1 + z2
if u_new <= max_u and z_new <= max_z:
y_new = tuple(y1[i] + y2[i] for i in range(c))
result[(u_new, z_new, y_new)] += c1 * c2
return dict(result)
def at_least_one_geometric(j, c, max_u, max_z):
"""
Compute u*y_j*z^j / (1 - u*y_j*z^j) truncated.
This represents "at least one card with value j".
"""
result = {}
for t in range(1, max_u + 1):
z_power = t * j
if z_power <= max_z:
y_tuple = tuple(t if idx == j - 1 else 0 for idx in range(c))
result[(t, z_power, y_tuple)] = 1
return result
def geometric_series(j, c, max_u, max_z):
"""
Compute 1 / (1 - u*y_j*z^j) truncated.
This represents "any number (including 0) of cards with value j".
"""
result = {}
for t in range(0, max_u + 1):
z_power = t * j
if z_power <= max_z:
y_tuple = tuple(t if idx == j - 1 else 0 for idx in range(c))
result[(t, z_power, y_tuple)] = 1
return result
def compute_theta(c, k, max_z=13):
"""
Compute Theta_k(c) using generating functions.
Returns: list of multisets (as sorted tuples) where A_k = c and sum <= max_z.
"""
max_u = k
# Start with: at least one card with value c
poly = at_least_one_geometric(c, c, max_u, max_z)
# Multiply by 1/(1-u*y_j*z^j) for each j < c
for j in range(1, c):
factor = geometric_series(j, c, max_u, max_z)
poly = multiply_poly(poly, factor, c, max_u, max_z)
# Extract [u^k][z^{<=max_z}] and decode multisets from y-tuples
multisets = []
for (u_pow, z_pow, y_tuple), coef in poly.items():
if u_pow == k and z_pow <= max_z:
elements = []
for j in range(c):
elements.extend([j + 1] * y_tuple[j])
elements.sort()
multisets.append(tuple(elements))
return sorted(multisets)
def labeling_coefficient(multiset, n):
"""
Compute Lab_n = n! / prod(m_j!)
"""
from collections import Counter
counter = Counter(multiset)
result = factorial(n)
for m in counter.values():
result //= factorial(m)
return result
def compute_K_ge_k(k, n_cards=4, max_z=13):
"""
Compute #(K >= k) for n_cards total cards.
"""
total_count = 0
details_by_c = {}
remaining_cards = n_cards - k
for c in range(1, 14):
theta_multisets = compute_theta(c, k, max_z)
if not theta_multisets:
continue
c_count = 0
for first_k_multiset in theta_multisets:
for extra in combinations_with_replacement(range(c, 14), remaining_cards):
full_multiset = tuple(sorted(first_k_multiset + extra))
lab = labeling_coefficient(full_multiset, n_cards)
c_count += lab
details_by_c[c] = {
'theta_multisets': theta_multisets,
'count': c_count
}
total_count += c_count
return total_count, details_by_c
# =============================================================================
# Main computation
# =============================================================================
print("=" * 70)
print("Computing #(K >= k) using Generating Functions")
print("=" * 70)
for k in [2, 3]:
print(f"\n{'=' * 70}")
print(f"K >= {k} (A_1 + A_2 + ... + A_{k} <= 13)")
print(f"{'=' * 70}")
result, details = compute_K_ge_k(k)
print(f"\n[Breakdown by A_{k} = c]")
print(f"{'c':<4} {'|Theta|':<8} {'Contribution':<14} {'Cumulative':<12} Theta_{k}(c)")
print("-" * 70)
cumsum = 0
for c in sorted(details.keys()):
cumsum += details[c]['count']
theta_size = len(details[c]['theta_multisets'])
multisets_str = str(details[c]['theta_multisets'])
print(f"{c:<4} {theta_size:<8} {details[c]['count']:<14} {cumsum:<12} {multisets_str}")
print(f"\n#(K >= {k}) = {result}")
print(f"P(K >= {k}) = {result}/28561 = {result/28561:.6f}")
# =============================================================================
# Summary
# =============================================================================
print("\n" + "=" * 70)
print("Summary")
print("=" * 70)
k2, _ = compute_K_ge_k(2)
k3, _ = compute_K_ge_k(3)
k4, _ = compute_K_ge_k(4)
print(f"\n{'k':<5} {'#(K >= k)':<12} {'P(K >= k)':<12} {'#(K = k)':<12} {'P(K = k)':<12}")
print("-" * 55)
total = 13**4
print(f"{1:<5} {total:<12} {total/total:<12.6f} {total - k2:<12} {(total-k2)/total:<12.6f}")
print(f"{2:<5} {k2:<12} {k2/total:<12.6f} {k2 - k3:<12} {(k2-k3)/total:<12.6f}")
print(f"{3:<5} {k3:<12} {k3/total:<12.6f} {k3 - k4:<12} {(k3-k4)/total:<12.6f}")
print(f"{4:<5} {k4:<12} {k4/total:<12.6f} {k4:<12} {k4/total:<12.6f}")
计算结果：
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
======================================================================
Computing #(K >= k) using Generating Functions
======================================================================
======================================================================
K >= 2 (A_1 + A_2 + ... + A_2 <= 13)
======================================================================
[Breakdown by A_2 = c]
c |Theta| Contribution Cumulative Theta_2(c)
----------------------------------------------------------------------
1 1 913 913 [(1, 1)]
2 2 2359 3272 [(1, 2), (2, 2)]
3 3 3289 6561 [(1, 3), (2, 3), (3, 3)]
4 4 3775 10336 [(1, 4), (2, 4), (3, 4), (4, 4)]
5 5 3889 14225 [(1, 5), (2, 5), (3, 5), (4, 5), (5, 5)]
6 6 3703 17928 [(1, 6), (2, 6), (3, 6), (4, 6), (5, 6), (6, 6)]
7 6 3048 20976 [(1, 7), (2, 7), (3, 7), (4, 7), (5, 7), (6, 7)]
8 5 1820 22796 [(1, 8), (2, 8), (3, 8), (4, 8), (5, 8)]
9 4 976 23772 [(1, 9), (2, 9), (3, 9), (4, 9)]
10 3 444 24216 [(1, 10), (2, 10), (3, 10)]
11 2 152 24368 [(1, 11), (2, 11)]
12 1 28 24396 [(1, 12)]
#(K >= 2) = 24396
P(K >= 2) = 24396/28561 = 0.854172
======================================================================
K >= 3 (A_1 + A_2 + ... + A_3 <= 13)
======================================================================
[Breakdown by A_3 = c]
c |Theta| Contribution Cumulative Theta_3(c)
----------------------------------------------------------------------
1 1 49 49 [(1, 1, 1)]
2 3 319 368 [(1, 1, 2), (1, 2, 2), (2, 2, 2)]
3 6 793 1161 [(1, 1, 3), (1, 2, 3), (1, 3, 3), (2, 2, 3), (2, 3, 3), (3, 3, 3)]
4 10 1399 2560 [(1, 1, 4), (1, 2, 4), (1, 3, 4), (1, 4, 4), (2, 2, 4), (2, 3, 4), (2, 4, 4), (3, 3, 4), (3, 4, 4), (4, 4, 4)]
5 13 1932 4492 [(1, 1, 5), (1, 2, 5), (1, 3, 5), (1, 4, 5), (1, 5, 5), (2, 2, 5), (2, 3, 5), (2, 4, 5), (2, 5, 5), (3, 3, 5), (3, 4, 5), (3, 5, 5), (4, 4, 5)]
6 12 1798 6290 [(1, 1, 6), (1, 2, 6), (1, 3, 6), (1, 4, 6), (1, 5, 6), (1, 6, 6), (2, 2, 6), (2, 3, 6), (2, 4, 6), (2, 5, 6), (3, 3, 6), (3, 4, 6)]
7 9 1170 7460 [(1, 1, 7), (1, 2, 7), (1, 3, 7), (1, 4, 7), (1, 5, 7), (2, 2, 7), (2, 3, 7), (2, 4, 7), (3, 3, 7)]
8 6 660 8120 [(1, 1, 8), (1, 2, 8), (1, 3, 8), (1, 4, 8), (2, 2, 8), (2, 3, 8)]
9 4 324 8444 [(1, 1, 9), (1, 2, 9), (1, 3, 9), (2, 2, 9)]
10 2 126 8570 [(1, 1, 10), (1, 2, 10)]
11 1 30 8600 [(1, 1, 11)]
#(K >= 3) = 8600
P(K >= 3) = 8600/28561 = 0.301110
======================================================================
Summary
======================================================================
k #(K >= k) P(K >= k) #(K = k) P(K = k)
-------------------------------------------------------
1 28561 1.000000 4165 0.145828
2 24396 0.854172 15796 0.553062
3 8600 0.301110 7885 0.276076
4 715 0.025034 715 0.025034
脚注
实际上，卡牌的点数并非均匀分布。点数【2】和【Q】的概率为$7/80$，而其它的点数概率为$3/40$。由于这种差异对计算结果的影响可以忽略不记，我们还是认为所有卡牌点数是均匀分布的。 ↩︎