我的 最初的 GSOC 提案 是实现修改 Mamba.jl 以使其能够拟合 Crosscat，这是一个通用贝叶斯模型，它使用嵌套在列方向 Dirichlet 聚类中的行方向 Dirichlet 聚类模型来拟合表格数据。该模型本身具有广泛的用途，但我选择该项目的主要原因是为了研究更通用的内容：改进对具有离散和连续参数混合的模型进行 MCMC 的工具。
最终，我未能完成完整的原始计划。但是，我确实在 Mamba 中实现了一个简单的 Dirichlet 1D 高斯混合模型。虽然该模型本身非常基础，但它确实需要成功地重构 Mamba 以启用可变数量的参数——这是一项大量的工作。基于这项工作，涉及多个维度和/或 Dirichlet 过程的更复杂的 Dirichlet 混合模型几乎是微不足道的，甚至像 Crosscat（或其改进版本）这样更复杂的东西在将来也更容易实现。我估计，虽然我提供的演示模型的实用性仅占完整 Crosscat 的一小部分，但我完成的实际工作大约占实现它的 75%。
项目动机
贝叶斯 MCMC 是统计学（以及几乎所有科学）工具箱中一个强大的通用工具。它允许人们灵活地构建模型，以捕捉已知动力学和未知参数、可测量数据和不确定噪声之间的相互作用，从而从数据中提取意义。
这个想法的强大之处在于它的灵活性。如果你可以编写似然模型，你至少可以尝试进行 MCMC。当然，计算和收敛问题可能会使实际操作变得困难，但至少在理论上，这个想法足够简单，可以适用于自动化。
目前，这方面最突出的工具是 Stan。Stan 的采样器是 NUTS（“无 U 转弯采样器”），它依赖于 HMC（即自动微分）来有效地探索后验分布，即使这些后验分布位于高维参数空间中的中等维流形上——对于旧的、HMC 之前的采样器来说，如果没有巧妙的特定问题技巧来使各个维度准独立，这将是无法实现的。而且 NUTS 甚至不需要其许多 HMC 兄弟姐妹所需的任何手动调整。
然而，Stan 的基本设计意味着它具有一些不太可能解决的弱点。首先，它使用专有的语言进行模型定义，这带来了所有由此带来的限制和摩擦。最终用户几乎肯定不想深入研究 Stan 的 C++ 代码来添加功能。其次，由于 NUTS 采样器是在其基础上构建的，因此在混合离散和连续参数的模型上使用它可能仍然很困难。
Julia，特别是 Mamba.jl，提供了超越这些弱点的方法。虽然目前还没有像 Stan 那样成熟，但 Mamba.jl 确实具有使用灵活语法构建模型的基本功能，该语法反映了统计学家思考模型的方式。Mamba 还已经拥有许多采样器选项，包括 NUTS。模型使用通用的 Julia 代码而不是单一用途的语言以富有表现力和灵活的方式定义；并且可以组合各种采样器，以便模型可以包含离散和连续参数。
但在本项目之前，Mamba 仅限于具有固定数量参数的模型。这为许多有用的模型类型关闭了大门。例如，在 Dirichlet 混合模型和其他聚类模型中，参数的数量取决于拟合模型在数据中找到的潜在聚类数量。这就是我的项目旨在填补的空白。
障碍
在执行此项目时，我遇到了几个意想不到的障碍。首先，有 这个（警告，血）。那是我在骨折后在急诊室做鬼脸的照片；我在医院里住了 3 天，又过了一周才停止服用止痛药并能够再次打字。总而言之，那次事故（经典的傻瓜在我骑车经过时打开车门，外加一场暴雨）可能让我损失了 2 周的工作时间。
此外，重构 Mamba 证明比我预期的要困难。我的计划是向许多基本的 Mamba 类型添加参数，以便能够在将参数存储在现有的固定大小数组结构中或在我的新设计的灵活大小结构中进行切换。趁此机会，我还向类型参数添加了类型参数，以放松对 Float64 模型参数的硬编码依赖关系，以便能够为 HMC 使用自动微分数字。对于我当时对 Julia 本身和 Mamba 包的专业知识水平来说，这相当高级；我收到了很多错误消息才真正理解了一些东西。（当然，现在我理解了，这似乎微不足道；但这确实是一场斗争，因为问题当然没有像我在下面展示的那样清晰地出现。）
例如
事实证明，即使你真的希望它们是协变的，Julia 类型也不是协变的。例如，即使
VecDictVariateVal和SymDictVariateVal是我的通用抽象类型DictVariateVal的微不足道的不同子类型，VecDictVariateVal{Float64} <: DictVariateVal{Float64}也不成立。这尤其令人困惑（至少对我这个相对新手来说），因为使用where子句，UnionAll类型可以是协变的。你可以编写
SomeType{<:Real}或SomeType{T} where T<:Real，但永远不能编写SomeType{T<:Real} where T；最后那个只是SomeType{false}，因为它是一个类别错误；类型变量本身永远不是实数的子类型。在这个简化的最小示例中，这似乎很明显，但相信我，有些情况很难看出来。在接下来的几天里，我将发布一些 Julia 问题（YouTube 链接，抱歉），并提供一些关于如何改进语法本身以及错误消息/警告（当你出错时）的建议。Mamba 的控制流有点难以理解。我发现探索像这样的大型现有包的一个好技巧是，在一个工作示例上运行 图形化配置文件浏览器；这将为你提供一个关于哪些函数调用哪些函数的有用信息。
总而言之，我花了 6 周多的时间才完成这个重构，而我最初乐观地计划在第一个月花费不到 1/3 的时间来完成它。
那么，模型在哪里？
我的工作在 https://github.com/jamesonquinn/Mamba.jl；主要在 这个分支。除了整体重构之外，关键文件还包括 Dirichlet 过程分布；可逆跳跃采样器（主要基于 Richardson 和 Green 1997，但也有一些简化以及一些更改，以便将其基于 Dirichlet 过程而不是每个聚类数量的权重的单独 Dirichlet 分布）；以及 演示示例模型。
这是我在该示例中使用的模型
for
**T ~ DP(α)**，一个维度为**N**的整数向量（聚类索引）。
**μt 独立同分布 ~ N(μ0, τ2)** 对于任意 **t**； **1/σt2 独立同分布 ~ Γ(α, β)** 对于任意 **t**； **β 独立同分布 ~ Γ(θ, ϕ)** 对于任意 **t**。
**α = 0.1**， **μ0 = ȳ**， **τ = 2sy**， **θ = ϕ = 0.1**。
以下是两个聚类中45个模拟数据点的结果，标准差为3，均值为±5。
如您所见，两条链大部分时间都至少有一个聚类分别围绕“正确”值，但偶尔也会出错。
进一步工作
这段代码虽然适用于示例模型，但尚未准备好合并到Mamba的主分支中。我最终没有时间完成几个清理步骤。
我已更新了“slice”和“slicesimplex”采样器以与新的数据结构一起使用。但是，在我工作中未使用到的其他采样器在gsocMNVP分支中目前已损坏；它们仍然尝试使用旧的数据结构。沿着与slice和slicesimplex采样器相同的思路更新它们，将是一个或多或少例行的任务——每个采样器大约需要一个小时的工作。
除了上面显示的轨迹图之外，诊断和绘图也没有更新。修复此问题是一个不太简单的任务，因为由于“标签问题”，大多数诊断都需要以某种方式重新考虑才能应用于狄利克雷模型。
一旦完成清理工作（几天工作）并完成合并，按照原始计划实现完整的Crosscat模型应该不会太难。乐观地估计，我觉得需要1-2周……这意味着现实地讲，可能需要4-6周更现实。无论如何，我实现的新数据结构将使这项工作主要成为实现统计算法的问题；数据和模型基础设施都已就绪。
通过结合NUTS和离散功能，我相信Mamba将开始在某些任务上真正优于Stan。它在赶上Stan的成熟度方面还有很长的路要走，但在解决“两种语言问题”方面，它给了我以及其他人继续这项工作的强烈动力。
致谢
我要感谢我的GSoC导师Benjamin Deonovic在他帮助我完成这个困难但有趣的项目中给予的帮助和理解。