用于计算分子生物学的 Python 工具
用于计算分子生物学的 Python 工具
Biopython 是生物信息学和计算生物学中一个非常受欢迎的库。Mocapy++ 是一个用于动态贝叶斯网络 (DBN) 中参数学习和推理的机器学习工具包PaluszewskiHamelryck2010，它对域中随机变量之间的概率关系进行编码。Mocapy++ 在 GNU 通用公共许可证 (GPL) 下免费提供，可从 SourceForge 获取。该库支持各种 DBN 架构和概率分布，包括来自方向统计的分布。值得注意的是，球体上的肯特分布和环面上的二元冯米塞斯分布，已被证明可用于构建蛋白质和 RNA 结构的概率模型。
这样一个非常有用且强大的库，在 TorusDBNBMTFKH2008、BasiliskHBPFJH2010、FB5HMMHKK2006PaluszewskiWinter2008 等项目中被成功使用，是长期努力的结果。最初的 Mocapy 实现可以追溯到 2004 年，从那时起，该库就被重写成 C++。然而，C++ 是一种静态类型化的编译编程语言，不利于快速原型设计。因此，目前 Mocapy++ 还没有为动态加载自定义节点类型提供任何规定，并且一种不需要修改和重新编译库就能插入新节点类型的机制也值得关注。这种插件接口将通过允许快速实现和测试新的概率分布来帮助快速原型设计，这反过来可以大幅减少开发时间和精力；用户将能够在不修改和重新编译的情况下扩展 Mocapy++。认识到这种需求，该项目（在此称为 MocapyEXT），旨在改进当前的 Mocapy++ 节点类型扩展机制，由 T. Hamelryck 提出。
MocapyEXT 项目主要是一项工程努力，旨在为 Mocapy++ 提供一个透明的 Python 插件接口，在该接口中，内置和动态加载的节点类型可以以统一的方式使用。此外，外部实现和动态加载的节点可以由用户修改，而这些更改不会导致客户端程序或伴随的 Mocapy++ 库重新编译。这将有助于快速原型设计，简化对现有代码的改编，并提高软件互操作性，同时对现有的 Mocapy++ 接口进行最小更改，从而促进对 MocapyEXT 引入的更改的顺利接受。
Justinas V. Daugmaudis ([email protected])
导师
Thomas Hamelryck
Eric Talevich
’'’了解 S-EM 和方向统计 ‘’’
’'’研究 Mocapy++ 用例 ‘’’
’'’研究 Mocapy++ 内部结构和代码 ‘’’
’'’设计 Mocapy++ 插件接口 ‘’’
’'’实现 Mocapy++ 插件模块 ‘’’
’'’使用模块化 Mocapy++ 架构进行实验 ‘’’
’'’第 1-2 周 [5 月 23 日 – 6 月 5 日] ‘’’
接口策略设计：插件 API 对 Mocapy++ 用户的“自然”感受至关重要，这一点再怎么强调也不为过。因此，大量时间将用于接口设计。
’'’第 3-5 周 [6 月 6 日 – 6 月 19 日] ‘’’
实现插件模块。
’'’第 6-7 周 [6 月 20 日 – 6 月 30 日] ‘’’
实现一些示例，展示 Mocapy++ 插件系统在实际中的应用；例如，如何使用外部实现的逻辑斯蒂分布。
’'’第 7-8 周（中期）[7 月 1 日 – 7 月 10 日] ‘’’
’'’第 9-10 周 [7 月 11 日 – 7 月 24 日] ‘’’
示例应用程序。应重点关注插件模块的类型反射功能。
’'’第 11-12 周 [7 月 25 日 – 8 月 7 日] ‘’’
更新文档以反映新功能。此外，记录示例，对代码进行任何清理工作等。计划的“收尾”日期。
’'’第 13-14 周 [8 月 8 日 – 8 月 21 日] ‘’’
向更广泛的受众介绍绑定，收集社区的意见和评论。
’'’第 15 周 [8 月 22 日] ‘’’
项目结束。
托管在 gSoC11 Mocapy 分支
对于不熟悉 Mocapy++ 术语的人来说，有必要了解 ESS 计算机和密度是什么，因为这些术语将在本文中使用。
在下表中，X 表示 ESS 计算机类，u 是X 的可变值。
| 表达式 | 返回类型 | 先决条件/后置条件 |
|---|---|---|
| u.construct(parent_sizes) | void | 定义了 ESS 的适当形状 |
| u.estimate(ess) | void | 将样本点添加到 ESS |
类 mocapy::BippoESS 是 ESS 计算机的示例模型。
在下表中，X 表示 Densities 计算机类，v 是X 的常量值，u 是X 的可变值。请注意，ptv 代表“父节点和当前节点的值”：父节点的值以及该方法所属节点的值。
| 表达式 | 返回类型 | 先决条件/后置条件 | |
|---|---|---|---|
| u.construct(parent_sizes) | void | 参数（均值、协方差、CPD 等）已初始化 | |
| u.estimate(ess) | void | 根据给定的 ESS 估计节点的参数 | |
| u.sample(ptv) | std::vector |
根据指示的父节点值返回一个样本。请注意，后续对 sample 成员函数的调用可能会返回不同的值，因此该操作是可变的 | |
| v.get_lik(ptv, log) | double | 返回似然度，P(子节点 | 父节点) |
| v.get_parameters() | std::vector<MDArray |
返回分布参数 | |
| v.get_node_size() | unsigned int | 返回节点大小 | |
| v.get_output_size() | unsigned int | 返回节点的输出大小 |
类 mocapy::BippoDensities 是 Densities 计算机的示例模型。
为了使客户端程序能够执行 Python 代码，有必要初始化脚本环境。这是通过调用 Py_Initialize() 函数完成的。解释器通过调用 Py_Finalize() 函数释放。
然而，重要的是要注意，Py_Initialize() 和 Py_Finalize() 调用之间的任何语句都可能会抛出异常。如果抛出异常，则必须在 try/catch 块中处理它，或者必须终止程序。考虑到这一点，前面的示例程序可以通过将 Py_Initialize() 和 Py_Finalize() 之间的语句封装在 try/catch 块中来变得更加健壮。Py_Finalize() 的安全性不容忽视。目前，Boost.Python 有一些全局（或函数静态）对象，它们的生存期阻止引用计数下降到零，直到 Boost.Python 共享对象被卸载。这会导致崩溃，因为当引用计数降到零时，没有解释器。这提出了一个问题，即这种初始化 Python 解释器的方法是否可以被认为是“易于使用”、“安全”甚至“非侵入性的”，而这些是 MocapyEXT 的主要设计原则。
解决此问题的方案既易于使用、安全又非侵入性，其优雅程度令人惊讶，并展示了 RAII（资源获取即初始化）习惯用法。
Python 解释器在进入主函数之前初始化，并在退出主函数之后释放。除了包含定义必要的 Python 插件包装器的头文件外，不需要额外的用户操作。
Python 解释器应链接到客户端程序（而不是 Mocapy++ 库）。MocapyEXT 的最终用户应负责为 Python 头文件和库提供额外的包含和/或库路径，这些路径是成功编译和链接客户端程序所必需的。
MocapyEXT 还会以一种方式实例化其静态数据成员，使其仅在库的编译单元中实例化一次。在 C++ 标准中，明确指出
“… 特别是，静态数据成员的初始化（以及任何相关的副作用）不会发生，除非静态数据成员本身以需要定义静态数据成员存在的方式使用。”[cpp_std2003, temp.inst]
用户不需要管理实例化的静态数据成员。
MocapyEXT 插件具有以下线程安全保证
对传递的参数容器进行并发可变访问需要同步，例如通过读写锁。可变访问包括更改容器中的值、调用使迭代器失效的成员函数、通过移动构造函数移动容器。
当导入包含 ESS 或密度定义的模块 foo 时，嵌入的 Python 解释器将在环境变量 PYTHONPATH 指定的目录列表和当前路径中搜索名为 foo.py 的文件。
在继续导入模块之前，会显式地将 sys.path 列表附加到当前路径。依赖于脚本位于客户端程序以外的其他目录中的用户应在 PYTHONPATH 环境变量中列出这些路径。
请注意，由于解释器还会在安装相关的默认路径中进行搜索，因此，包含 ESS 或 Densities 定义的文件的名称不能与标准模块的名称相同。
最小接口导入示例展示了三种插件类的用法：densities_plugin、ess_plugin 和 aggregate_plugin。densities_plugin 和 ess_plugin 的预期用法模式是 ESS 和 Densities 类型分别在不同的文件中实现，即每个文件包含一个类。每个类都在其独立的解释器环境中加载。aggregate_plugin 的目的是简化已加载类型的管理，并优化嵌入式 Python 解释器的资源分配，因为它只为 ESS 和 Densities 类型创建了一个嵌入式解释器实例。
基本上，用户提供 Python 库的名称、类的名称以及构建类模型所需的事实参数列表。MocapyExt 库则返回新类实例。
给定的参数列表（最多支持 N = 6 个参数）将参数化 Densities 对象的初始化。参数将转发到 Python API，并使用 Boost.Python 库提供的类型转换功能进行转换。该机制是通用的，用户可以将任何任意类型 T 转换为其在 Python 中的等效类型。
参数通过常量引用转发。此方法接受并转发任意类型的参数，但代价是始终将参数视为常量。此解决方案通常用于构造函数参数。这种方法的一个特殊问题是无法形成对函数类型的常量引用，但这个缺陷（在核心问题 295 中已得到解决）实际上对我们有利，因为它可以防止恶意地尝试将函数传递给 Python 解释器。
对 ess 和 dens 成员函数的任何进一步调用都会自动调用 Python 中相应的类方法。
ESS 计算器和 Densities 对象（分别为 ess 和 dens）的生命周期可能超过其各自工厂插件的生命周期。ess 和 dens 对象保留指向 Python 解释器的引用计数指针，因此有效 Python 解释器实例将一直存在，直到 ess 和 dens 对象被销毁。
这是一个简单的示例，展示了初始化插件节点的两种不同方法。在本例中，程序中使用的模块 plugin_tests 实现了一个固定长度为 1 的虚拟离散节点，并且始终返回 [0,] 作为采样值。
plugin_node_type 是 ChildNode 类模板的 typedef。这里我们还注意到节点的生命周期独立于插件工厂的生命周期。
插件节点是可流式的，即它可以通过运算符输出到 std::ostream
输出内容与repr(Z) 相同，其中Z 是 ESS 计算器或 Densities 对象。
MocapyEXT 保留了旧的 Mocapy++ ChildNode 类模板序列化行为，这意味着对于不需要序列化的 ESS 计算器（包括 kentess、gaussianess、poissoness 等），不会进行序列化。但是，MocapyEXT 必须序列化实现 ESS 计算器的 Python 插件，以便保留模块和类的名称，以便 ESS 计算器以后可以被反序列化。
还需要注意的是，反序列化后，ESS 计算器和 Densities 将不会共享相同的 Python 解释器实例，即使它们在同一个模块中实现。
衡量插件接口的成本与 ESS 和 Densities 计算器本身的实现相比如何，很有意思。
我们测量名为 N、S 和 A 的测试的相对性能。所有测试都调用 ESS 和 Densities 计算器的成员函数，尽管没有进行任何特定计算。
N 测试代表 N(ative)，这意味着测试中使用了 ESS 和 Densities 计算器的纯 C++ 实现。
S 测试代表 S(eparate)；Python 类在单独的 Python 解释器实例中加载。
A 测试代表 A(ggregate)；ESS 和 Densities 计算器在同一个 Python 解释器实例中加载。
Conf. int. 4.53e+04 5.03e+04 5.53e+04
基于 Wilcoxon 配对置信区间的加速百分比。
Func. A vs Func. B Minimum Median Maximum (% faster) (% faster) (% faster) N vs S 717 719 724 N vs A 718 725 730 A vs S -0.673 -0.251 0.0774
测试结果表明，使用 MocapyEXT 插件接口会带来一些性能损失。然而，必须注意的是，测试还执行了 ESS 和 Densities 实例的重复构建。
Conf. int. 1.94e+04 2.15e+04 2.37e+04
基于 Wilcoxon 配对置信区间的加速百分比。
Func. A vs Func. B Minimum Median Maximum (% faster) (% faster) (% faster) N vs S 240 242 243 N vs A 239 243 245 A vs S -1.24 -0.51 -0.25
重复测试，不构建 ESS 和 Densities 对象，结果表明，通过 MocapyEXT 接口调用方法的速度仅比调用本机实现的 ESS 和 Densities 对象的成员函数慢约 2.5 倍。很明显，在实际场景中，节点构建、采样、似然计算等操作更可能占运行时间的大部分，而不是方法本身的调用。基本上，这些测试表明，成员函数调用所花费的时间远小于像节点构建这样的“轻量级”操作。原则上，现在可以编写相当通用的算法。考虑以下示例