Bayesian Modeling of Collatz Stopping Times: A Probabilistic Machine Learning Perspective

本文从概率机器学习视角出发，通过贝叶斯分层负二项回归与基于奇数块分解的生成式近似两种模型，对 $n \le 10^7$ 范围内的 Collatz 总停止时间进行了建模分析，发现低阶模结构（如 $n \bmod 8$）是解释其异质性的关键因素，且统计回归模型在预测性能上优于生成式近似。

要点

这篇论文就像是在用**“天气预报”和“乐高积木”两种不同的视角，去研究一个数学界著名的谜题——“考拉兹猜想”**（也叫 $3n+1$ 问题）。

为了让你轻松理解，我们把这篇充满数学公式的论文，翻译成几个生动的故事：

1. 背景：那个让人抓狂的“数字游戏”

想象你有一个数字游戏：

• 如果数字是偶数，就把它除以 2。
• 如果数字是奇数，就把它乘以 3 再加 1。
• 一直重复，直到数字变成 1 为止。

考拉兹猜想说：无论你怎么开始（只要是正整数），最后一定能变成 1。虽然没人能证明这一点，但数学家们发现，如果你随机挑一个很大的数字玩这个游戏，它变成 1 需要的步数（论文里叫 $\tau(n)$）是有规律可循的。

这篇论文就是去研究：“如果我们随机挑一个数，它大概要走多少步才能回家（变成 1）？”

---

2. 数据发现：步数像“长尾怪兽”

作者计算了从 1 到 1000 万的所有数字的步数。他们发现：

• 大多数数字走的步数差不多，集中在某个平均值附近。
• 但是，极少数数字会走非常非常长的路（比如几千步），就像怪兽的长尾巴一样。
• 而且，步数的多少不仅取决于数字有多大，还取决于这个数字**“长什么样”**（比如它除以 8 余几）。

这就像是在统计**“人们去上班的时间”**：大部分人 30 分钟到，但偶尔有人因为堵车要 3 小时，而且这种“堵车”往往和住在哪个街区（模 8 的余数）有关。

---

3. 方法一：聪明的“天气预报员” (贝叶斯负二项回归)

作者首先请来了一个**“天气预报员”**（统计模型）。

• 它的逻辑：它不看数字内部怎么变化的，它只看两个简单的特征：
  1. 数字有多大（取对数后）：数字越大，通常步数越多，就像路程越远，开车时间越长。
  2. 数字的“身份证”（除以 8 的余数）：就像住在不同街区，通勤时间会有系统性差异。
• 它的超能力：它使用了一种叫**“负二项分布”的数学工具。这就像天气预报员知道，虽然平均时间是 30 分钟，但时间波动很大（有的 10 分钟，有的 60 分钟），所以它不会只报一个数，而是报一个“概率范围”**（比如：90% 的概率在 20-40 分钟之间）。
• 结果：这个“天气预报员”非常准！它能很好地预测出大多数数字的步数，并且能给出一个“不确定性”的区间。在预测测试中，它的表现是最好的。

比喻：这就像你问天气预报员“明天去上班要多久？”，他看着你的住址和距离，告诉你：“大概率 30 分钟，但考虑到路况波动，可能是 25 到 45 分钟。”

---

4. 方法二：拆解“乐高积木” (机械生成模型)

既然“天气预报员”很准但有点“黑箱”（不知道具体为什么），作者又造了一个**“乐高机械师”**（生成模型）。

• 它的逻辑：它把考拉兹游戏拆解成一个个**“积木块”**。
  • 当一个奇数变成 $3n+1$ 后，它一定会变成偶数，然后被连续除以 2 很多次，直到再次变成奇数。
  • 作者把这一连串“除以 2"的过程看作一个**“积木块”**（Block）。
  • 这个积木块有多长（被除以 2 几次），作者假设它是随机的。
• 它的进化：
  • 初版：假设积木块长度完全随机（像抛硬币）。结果发现，这太粗糙了，预测不准。
  • 升级版：作者发现，积木块的长度其实和数字**“除以 8 余几”**有关！就像不同街区的交通拥堵模式不同。于是，他们给积木块加上了“条件”：如果是余 1 的数，积木块长这样；如果是余 3 的数，积木块长那样。
• 结果：这个“机械师”虽然能解释为什么会有长尾巴（因为有时候积木块特别长），但在预测的准确度上，还是不如那个“天气预报员”。不过，它让我们明白了数字背后的机械原理。

比喻：这就像你不仅想知道“要多久”，还想知道“为什么”。机械师说：“因为你在 A 街区，红绿灯多，所以每次都要多等 3 个积木块的时间。”

---

5. 最终结论：谁赢了？

作者把两个模型放在测试集上比试（就像让两个司机在陌生城市比赛）：

1. 预测冠军：“天气预报员”（统计模型）。它最擅长猜中具体的步数，因为它利用了数据的统计规律，而且考虑了波动性。
2. 解释冠军：“乐高机械师”（生成模型）。虽然它猜得没那么准，但它揭示了**“模 8 余数”**这个关键因素。它告诉我们：数字的步数差异，很大程度上是因为它们属于不同的“模 8 家族”。

一句话总结：
这篇论文告诉我们，虽然考拉兹猜想的数学证明很难，但我们可以用统计学（像天气预报一样精准预测）和概率生成（像拆解积木一样理解机制）两种方法，把那个神秘的数字游戏变得有迹可循。特别是，**“数字除以 8 的余数”**是决定它要走多远的关键秘密。

这就好比，虽然我们无法预知每一辆车的精确到达时间，但我们可以通过统计规律和交通结构，非常准确地描述整个城市的交通状况。

技术摘要

论文技术总结

1. 研究问题与背景

• 核心对象：Collatz 猜想中的总停止时间 $\tau(n)$，即从整数 $n$ 开始，经过 Collatz 映射 $T(n)$（$n$ 为偶数除以 2，$n$ 为奇数则 $3n+1$）直到首次到达 1 所需的步数。
• 研究目标：作者并不试图证明 Collatz 猜想，而是将 $n$ 视为随机变量，研究在大规模数据（$N=10^7$）下 $\tau(n)$ 的经验分布规律及其算术异质性（Arithmetic Heterogeneity）。
• 数据特征：
  • $\tau(n)$ 是一个高度偏斜且过度离散（Overdispersed）的计数变量。
  • 存在显著的算术结构：$\tau(n)$ 的分布随 $n$ 的变化呈现带状结构，且方差随均值增加而增大（异方差性）。
  • 经验离散度比率 $R = \widehat{\text{Var}}(\tau) / \widehat{\mathbb{E}}[\tau] \approx 24.56 \gg 1$，排除了泊松分布作为噪声模型的可能性。

2. 方法论

作者提出了两种互补的建模方法：一种基于统计回归的预测模型，另一种基于机制生成的近似模型。

方法一：贝叶斯分层负二项回归 (Bayesian Hierarchical NB2-GLM)

• 模型选择：鉴于数据的过度离散特性，采用负二项分布 (Negative Binomial, NB2) 作为似然函数。NB2 允许方差随均值二次增长 ($\text{Var} = \mu + \alpha\mu^2$)。
• 协变量设计：
  • 尺度项：$\log n$。实验表明 $\tau(n)$ 随 $n$ 的增长近似于对数线性关系。
  • 算术结构项：$n \pmod 8$。利用模 8 的余数类作为分类特征，捕捉低阶模结构带来的异质性。
• 分层结构：
  • 构建广义线性模型 (GLM)：$\log \mu_n = \beta_0 + \beta_{\log} \log n + u_{n \pmod 8}$。
  • 引入随机截距 $u_r \sim \mathcal{N}(0, \sigma_u^2)$ 对应 8 个余数类。这种部分池化（Partial Pooling）策略既能捕捉各类别的特异性，又能防止过拟合。
• 推断：使用 PyMC 框架和 NUTS 采样器进行贝叶斯推断，通过弱信息先验估计后验分布。

方法二：基于奇数块分解的随机生成模型 (Mechanistic Odd-Block Generative Model)

• 机制基础：利用 Collatz 映射的奇数块分解。对于奇数 $m$，$3m+1$可分解为$2^{K(m)} m'$，其中$K(m) = v_2(3m+1)$是 2 的幂次。一次奇数更新后跟随$K(m)$ 次除以 2 的操作。
• 随机化近似：
  • 将确定性的 $K(m)$ 替换为随机变量序列 $K_j$，其概率质量函数 (PMF) 为 $p_k$。
  • 生成停止时间：$\tau_{\text{gen}}(n) = v_2(n) + \sum (1 + K_j)$。
• 校准与改进：
  • 基准假设：经典启发式认为 $K$ 服从几何分布 ($P(K=k) \approx 2^{-k}$)。
  • 全局校准 (G2)：使用狄利克雷先验 (Dirichlet Prior) 从数据中估计全局 $p_k$。
  • 条件校准 (G3)：将 $p_k$ 的分布条件化于 $m \pmod 8$，即针对不同余数类估计不同的块长度分布。
• 实现细节：使用狄利克雷 - 多项式共轭更新来估计参数，并引入奇数投影算子 $\lfloor \cdot \rceil_{\text{odd}}$ 处理随机更新中可能产生的偶数中间态。

3. 实验设置与评估

• 数据集：$N=10^7$，使用动态规划（记忆化搜索）和 Numba JIT 加速计算 $\tau(n)$。
• 划分：训练集 $N_{\text{fit}}=50,000$，测试集 $N_{\text{test}}=50,000$（均匀随机抽取）。
• 评估指标：
  • 对数预测分数 (Log Predictive Score)：基于严格评分规则，衡量模型对真实观测值的概率赋值能力（越高越好）。
  • 1-Wasserstein 距离 (W1)：衡量预测分布与经验分布的整体形状差异（越低越好）。

4. 主要结果

• 预测性能对比（见表 2）：
  • NB2-GLM (M3) 表现最佳：对数得分为 -272,911.95，W1 距离为 3.20。它显著优于生成模型，能更准确地预测停止时间的分布。
  • 全局生成模型 (G2)：表现最差（得分 -1,165,983），表明简单的几何分布假设不足以捕捉复杂的算术异质性。
  • 条件生成模型 (G3)：引入 $m \pmod 8$ 条件后，性能大幅提升（得分 -1,079,086，W1 降至 5.43），但仍不及回归模型。
• 分布拟合：
  • NB2-GLM 能很好地拟合主体分布，仅轻微高估了极右尾。
  • 条件生成模型 (G3) 通过后验预测检查 (PPC) 显示，条件化显著减少了均值偏差并改善了主体分布的吻合度，证明了低阶模结构是异质性的关键驱动因素。

5. 关键贡献

1. 统计视角的量化：首次从概率机器学习的角度，系统量化了 Collatz 停止时间的过度离散性和算术异质性，并证明了负二项回归的有效性。
2. 模块化结构的发现：通过分层回归和条件生成模型，明确证实了模 8 余数类（$n \pmod 8$）是解释 $\tau(n)$ 变异性的关键低维特征。
3. 机制与统计的桥梁：
   • 回归模型中的“随机截距”在生成模型中转化为显式的“条件概率分布”。
   • 展示了如何将确定性的数论问题转化为可校准的随机生成过程。
4. 方法论创新：提出了一种结合贝叶斯分层模型（用于预测和不确定性量化）与机制生成模型（用于解释物理/数论结构）的互补分析框架。

6. 意义与展望

• 理论意义：虽然未解决 Collatz 猜想，但该研究揭示了 Collatz 动力学中隐藏的统计规律，表明其停止时间并非完全随机，而是受低阶模结构强约束。
• 应用价值：为处理具有复杂算术结构的确定性动力系统提供了新的分析范式（Working Likelihood 视角）。
• 未来工作：
  • 将条件结构扩展至更高次幂的 2（如 $n \pmod{16}$）。
  • 在生成模型中引入显式的状态依赖性（State-dependence）。
  • 利用基于似然的校准方法，进一步对齐生成模型的评分规则表现与机制可解释性。

---

总结：这篇论文通过大规模数据实验和贝叶斯建模，成功地将 Collatz 停止时间问题转化为一个可预测的统计问题。它证明了简单的统计模型（NB2-GLM）在预测精度上优于复杂的机制生成模型，但后者在解释“为什么”存在异质性方面具有独特价值，特别是揭示了模 8 结构的核心作用。