用二元泊松模型预测2022世界杯8强

网上有很多文章用双泊松（Double Poisson）模型来预测世界杯比赛结果。但是双泊松模型有一个严重的缺陷，那就是它假设比赛中两队的比分是条件独立的。而我们都知道，在对抗性比赛中，两队的比分是存在关联的，因为两队都会根据场上的比分形势调整策略。比如足球比赛，当主队1:0领先，且距离比赛结束只剩10分钟时，落后的客队会孤注一掷，甘愿冒更大风险去争取平局。但如果主队3:0甚至4:0领先时，领先的主队可能会稍微放松下来，甚至教练会用新人换下主力，此时落后的客队更容易进1球（甚至主队会礼貌性让球）。所以比赛中两队比分是相关的，这种相关性可以通过依赖性参数来描述。

二元泊松（Bivariate Poisson）模型可以度量两队比分的依赖性参数，用二元泊松模型对比赛进行的预测准确率更高，在1/8决赛已经进行的4场比赛中，二元泊松模型预测正确率100%。

需要完整源代码的朋友可以关注私信我，或者评论留言索取。

二元泊松模型

考虑随机变量Xr,r∈{1,2,3}X_r, r \in \{1,2,3\}Xr​,r∈{1,2,3}服从独立泊松分布，其参数${\lambda }_r>0$，那么随机变量$X=X_1+X_3$和$Y=Y_2+Y_3$服从二元泊松分布$BP({\lambda }_1,{\lambda }_2,{\lambda }_3)$，其联合概率为：

PX,Y(x,y)=Pr(X=x,Y=y)=exp⁡{−(λ1+λ2+λ3)}λ1xx!λ2yy!×∑k=0min(x,y)(xk)(yk)k!(λ3λ1λ2)k\begin{aligned} P_{X,Y}(x,y)&=Pr(X=x,Y=y)\\ &=\exp\{-(\lambda_1+\lambda_2+\lambda_3)\}\frac{\lambda_1^x}{x!}\frac{\lambda_2^y}{y!}\times\sum_{k=0}^{min(x, y)}\begin{pmatrix}x\\k\end{pmatrix}\begin{pmatrix}y\\k\end{pmatrix}k!(\frac{\lambda_3}{\lambda_1\lambda_2})^k \end{aligned} PX,Y​(x,y)​=Pr(X=x,Y=y)=exp{−(λ1​+λ2​+λ3​)}x!λ1x​​y!λ2y​​×k=0∑min(x,y)​(xk​)(yk​)k!(λ1​λ2​λ3​​)k​

边缘概率上，每个随机变量服从泊松分布，期望$E(X)={\lambda }_1+{\lambda }_3$，$E(Y)={\lambda }_2+{\lambda }_3$，协方差$cov(X,Y)={\lambda }_3$。这里${\lambda }_3$就是度量两队进球之间相关性的指标。如果${\lambda }_3=0$，则两个变量是条件独立的，此时二元泊松分布退化成两个独立泊松分布的乘积，即双泊松模型。

我们用$(x_n,y_n)$表示第$n$场比赛主队和客队的比分，那么带比分相关性的二元泊松模型的一般形式为：

$$\begin{array}{rl}{X}_n,Y_n\mid {\lambda }_{1n},{\lambda }_{2n},{\lambda }_{3n}& \sim \text{BivPoisson}({\lambda }_{1n},{\lambda }_{2n},{\lambda }_{3n})\\ \log ({\lambda }_{1n})& =\mu +home+at{t}_{h_n}+de{f}_{a_n}\\ \log ({\lambda }_{2n})& =\mu +at{t}_{a_n}+de{f}_{h_n}\\ \log ({\lambda }_{3n})& ={\beta }_0+{\gamma }_1{\beta }_{h_n}^{home}+{\gamma }_2{\beta }_{a_n}^{away}+{\gamma }_3\beta w_n\end{array}$$

其中:

${\lambda }_{1n},{\lambda }_{2n}$分别表示主队和客队的进球率；

$\mu$表示截距项；

$home$代表主场效应，众所周知在自己的主场踢球是由很多优势的（更熟悉场地、更适应气候、更热情的球迷…）;

$at{t}_t$ 和 $de{f}_t$ 分别表示每只球队的进攻能力和防守能力，公式中下标$t$写作$h_n,a_n$，分别表示主队和客队的第$n$场比赛；

${\beta }_{h_n}^{home}$ 和 ${\beta }_{a_n}^{away}$ 分别代表主队和客队的参数；

$w_n$是第$n$场比赛的协方差向量，用于建模协方差项；

$\beta$ 是对应的回归系数向量；

参数${\gamma }_1,{\gamma }_2,{\gamma }_3$是取值为0或1的二进制指示符，可激活线性预测的不同来源。因此，当${\gamma }_1={\gamma }_2={\gamma }_3=0$时，我们得到常数协方差；而当$({\gamma }_1,{\gamma }_2,{\gamma }_3)=(1,1,0)$时，我们假设协方差仅取决于团队的参数，跟比赛协方差无关。

为了实现模型的可识别性，攻击/防御参数都是零和的：

$$\sum _{t=1}^{T}at{t}_t=0,\sum _{t=1}^{T}de{f}_t=0$$

这条约束大量出现在足球文献中，即假设第$T$只球队的能力等于其他球队能力和的负值，即能力的零和性。

$$at{t}_T=-\sum _{t=1}^{T-1}at{t}_t,de{f}_T=-\sum _{t=1}^{T-1}de{f}_t$$

有了上面的理论基础，我们就可以构造并训练我们的模型了。

最大似然法

给定参数向量 θ=({attt,deft,t=1,…,T},μ,home,βhnhome,βanaway,β0,β)\theta = (\{att_t, def_t, t=1,\dots, T\}, \mu, home, \beta_{h_n}^{home}, \beta_{a_n}^{away}, \beta_0, \beta)θ=({attt​,deft​,t=1,…,T},μ,home,βhn​home​,βan​away​,β0​,β) ，二元泊松模型的似然函数为：

L(θ)=∏n=1Nexp⁡{−(λ1n+λ2n+λ3n)}λ1nxnxn!λ2nynyn!×∑k=0min(xn,yn)(xnk)(ynk)k!(λ3nλ1nλ2n)kL(\theta) = \prod_{n=1}^N\exp\{-(\lambda_{1n}+\lambda_{2n}+\lambda_{3n})\}\frac{\lambda_{1n}^{x_n}}{x_n!}\frac{\lambda_{2n}^{y_n}}{y_n!}\times\sum_{k=0}^{min(x_n, y_n)}\begin{pmatrix}x_n\\k\end{pmatrix}\begin{pmatrix}y_n\\k\end{pmatrix}k!(\frac{\lambda_{3n}}{\lambda_{1n}\lambda_{2n}})^k L(θ)=n=1∏N​exp{−(λ1n​+λ2n​+λ3n​)}xn​!λ1nxn​​​yn​!λ2nyn​​​×k=0∑min(xn​,yn​)​(xn​k​)(yn​k​)k!(λ1n​λ2n​λ3n​​)k

最大似然参数估计只要找到最大似然估计函数 $\hat{\theta }$ 即可：

$$\hat{\theta }=\underset{\theta \in \Theta }{\mathrm{argmax}}L(\theta )$$

这个过程可以用求导解决：

$$l^{\prime }(\theta )=0$$

Wald检验+偏置信区间也可以构造最大似然估计$\hat{\theta }$，95%的Wald类型区间满足：

$$\hat{\theta }\pm 1.96se(\hat{\theta })$$

最大似然法仅适用于静态模型，因为静态模型复杂度不高。当使用动态模型时，随着参数空间的增长，最大似然法的计算量会非常大，且输出结果也不可靠。所以我们会更关注另一种方法——贝叶斯法。

贝叶斯法

贝叶斯分析的目标是从联合后验分布$\pi (\theta \mid D)$中得出推断结论，这里$D=(x_n,y_n{)}_{n=1,\dots ,N}$ 表示$N$场比赛的观测数据集合。联合后验分布满足：

$$\pi (\theta \mid D)=\frac{p(\theta \mid D)\pi (\theta )}{p(D)}\propto p(D\mid \theta )\pi (\theta )$$

其中 $p(D\mid \theta )$ 是模型采样分布（与似然函数成比例），$\pi (\theta )$ 是 $\theta$ 的联合先验分布。$p(D)={\int }_{\Theta }p(D\mid \theta )\pi (\theta )d\theta$ 是不依赖于$\theta$的边际似然。

在大多数情况下，$\pi (\theta \mid D)$不具有闭合形式，因此，我们需要通过模拟对其进行近似。处理这种情况的最主流的方法是马尔可夫链蒙特卡罗模拟。

就推断结论而言，我们通常对单一参数的边际后验分布（后验均值、中值、可信区间等）感兴趣。我们可以将上述二元泊松模型的后验分布公式写成：

$$\pi (\theta \mid D)\propto \pi (\theta )\prod _{n=1}^N\text{BivPoisson}({\lambda }_{1n},{\lambda }_{2n},{\lambda }_{3n})$$

其中$\pi (\theta )=\pi (att)\pi (def)\pi (\mu )\pi (home)\pi ({\beta }_{h_n}^{home})\pi ({\beta }_{a_n}^{away})\pi ({\beta }_0)\pi (\beta )$ 是先验独立参数分量假设下的联合最优分布。

标准方法是将一些信息量较弱的先验分布分配给团队的某些特定能力。这些参数在两个常见（先验）分布中可以交换：

$$\begin{array}{rl}& {\text{att}}_t\sim N({\mu }_{att},{\sigma }_{att})\\ & {\text{def}}_t\sim N({\mu }_{def},{\sigma }_{def}),t=1,\dots ,T\end{array}$$

这里面 $m{u}_{att},{\sigma }_{att},{\mu }_{def},{\sigma }_{def}$ 都是超参，通过向其余参数分配一些信息量较弱的先验来完成模型公式化。

用模型预测8强

我用2018-2022年期间的所有国际比赛，以及2022年世界杯小组赛的数据，通过Hamilton Monte Carlo采样、2000次迭代训练了一个模型。其思路是提供一个动态预测场景：在每个比赛日结束时，重新调整模型以预测剩余的比赛。对2022年世界杯淘汰赛第一轮的16场比赛的预测，每只球队的动态先验会集中关注参赛球队小组赛3场比赛的进攻防守数据，前三场比赛将作为一个单独独立的时间段，而不是作为三个不同的时间段来考虑。这里有一个因素需要考虑，就是小组赛的最后一场比赛中，一些球队已经晋级，他们可能没有派主力阵容或没有全力去踢，所以对这样的球队（比如葡萄牙），在先验数据上会进行一些微调。

用上述模型的后验预测分布的后验匹配概率见下表：

强队	弱队	胜	平	负	最可能结果
荷兰	美国	0.471	0.298	0.231	1-0 (0.164)
阿根廷	澳大利亚	0.665	0.246	0.089	1-0 (0.208)
法国	波兰	0.620	0.235	0.145	1-0 (0.145)
英格兰	塞内加尔	0.632	0.235	0.132	1-0 (0.162)
日本	克罗地亚	0.363	0.288	0.349	1-1 (0.122)
巴西	韩国	0.740	0.190	0.070	1-0 (0.171)
西班牙	摩洛哥	0.562	0.270	0.168	1-0 (0.173)
葡萄牙	瑞士	0.486	0.273	0.241	1-0 (0.141)

上表中，胜/负针对的都是强队，最可能结果后面的括号表示后验概率。

更加直观可视化的结果见下图：

上图中越深的颜色表示越可能出现的结果，x轴是强队，y轴是弱队。

从预测结果看，目前已经结束的4场1/8决赛全部预测正确:

比赛	预测	实际结果
荷兰-美国	✅荷兰胜 (0.471)	3:1 荷兰胜
阿根廷-澳大利亚	✅阿根廷胜 (0.665)	2:1 阿根廷胜
法国-波兰	✅法国胜 (0.620)	3:1 法国胜
英格兰-塞内加尔	✅英格兰胜 (0.632)	3:0 英格兰胜

---

12月6日更新

昨晚今晨进行的日本-克罗地亚和巴西-韩国的两场八分之决赛中，日本-克罗地亚120分钟踢平，点球大战中日本队2:4负于克罗地亚。这个结果与模型的预测相同。模型预测120分钟的最可能结果是1-1战平。另一场巴西-韩国毫无悬念，巴西4:1战胜韩国，这与模型预测巴西有74%的概率获胜一致。

比赛	预测	实际结果
日本-克罗地亚	✅120分钟战平 日本胜(0.363) - 平局(0.288 ) - 克罗地亚胜(0.349) 最可能结果：1-1(0.122)	2:4 克罗地亚胜（点球）
巴西-韩国	✅巴西胜 (0.740)	4:1 巴西胜

---

需要完整源代码的朋友可以关注私信我，或者评论留言索取。