企业官网建设流程全解析

• 科学计算
• 科研
• 机器学习

导读

rdkit.ML.Cluster.Clustering是 RDKit 中由 C++/Fortran 实现并导出给 Python 的层次聚类核心扩展模块，负责执行 Ward、单链（SLINK）、全链（CLINK）、组平均（UPGMA）等经典凝聚式层次聚类算法。本文以该模块为骨架，结合仓库中 Fortran 算法源码、C++ 封装、Python 接口与配套工具链，完整讲解其调用方式、距离矩阵约定、聚类树数据结构，以及如何基于分子指纹完成化合物集的实际聚类任务，帮助读者从"会调 API"深入到"懂原理、可排错、能扩展"。

一、模块定位：文档与源码对应关系

Docs/Book/source/rdkit.ML.Cluster.Clustering.rst是一份 Sphinxautomodule文档页，它声明了对rdkit.ML.Cluster.Clustering模块的自动文档化（:members:、:undoc-members:、:show-inheritance:三个指令分别表示生成公开成员、未文档化成员与继承关系的文档）。

该模块在仓库中的真实身份是一个编译型 Python 扩展模块：

• C++ 封装源码：Code/ML/Cluster/Murtagh/Clustering.cpp（Boost.Python 版本）与 Code/ML/Cluster/Murtagh/nbClustering.cpp（nanobind 版本）；
• 底层算法核心：Fortran 编写的 Code/ML/Cluster/Murtagh/hc.f（凝聚过程）与 Code/ML/Cluster/Murtagh/hcdriver.f（数据驱动层）；
• Python 高层封装：rdkit/ML/Cluster/Murtagh.py。

从模块依赖关系看，Murtagh.py通过from rdkit.ML.Cluster.Clustering import MurtaghCluster, MurtaghDistCluster导入该扩展，二者合起来构成了 RDKit 的 Murtagh 层次聚类接口。因此本文讨论的"Clustering 模块"实际指代这条完整的技术链路：Fortran 算法 → C++ 桥接 → Python 接口 → 聚类树对象。

二、算法内核：Murtagh 凝聚式层次聚类（hc.f）

hc.f文件头部注释明确说明了算法来源：F. Murtagh，ESA/ESO/STECF，Garching，1986 年 2 月。它实现的是经典的凝聚式（agglomerative）层次聚类：从 N 个单点簇开始，每轮合并最相似的两个簇，重复 N-1 次直到只剩一个簇，同时记录完整的合并历史。

2.1 核心数据结构

HC子程序的参数设计体现了算法的关键约定：

相异度矩阵的索引映射由IOFFSET函数完成：将对称矩阵的上三角（或下三角）元素压平为一维向量，即i<j时ioffset = (j-1)*(j-2)/2 + i。这一压缩存储格式贯穿整个模块，Python 侧_LookupDist与Resemblance.py均遵循同一约定（详见第五节）。

2.2 七种凝聚准则（IOPT）

hc.f中按IOPT分支实现了七种经典的簇间距离（链接）准则，对应 rdkit/ML/Cluster/Murtagh.py 中定义的常量：

其中m(·)为簇成员数（MEMBR），d(·,·)为相异度。值得注意的实现细节：

• Ward 法（IOPT=1）在初始化阶段会将所有距离除以 2，源码注释说明这是因为该方法的合并准则以方差而非距离定义；
• 单链、全链、组平均、McQuitty、中位数、质心法等公式在每次凝聚后就地更新被合并簇与其他簇之间的距离（DISS(IND1)覆盖写），实现 O(N²) 空间与接近 O(N³) 的时间复杂度。

2.3 算法流程

HC子程序按如下循环推进（对应源码中的标号400 CONTINUE主循环）：

1. 初始化：MEMBR(i)=1，FLAG(i)=.TRUE.，NCL=N；
2. 建立最近邻表NN/DISNN：对每个 i，扫描 j>i 找到最小相异度；
3. 全局扫描NN表，找出当前最小相异度DMIN对应的簇对(IM, JM)；
4. 记录合并：IA(N-NCL)=I2、IB(N-NCL)=J2、CRIT(N-NCL)=DMIN，其中I2=min(IM,JM)、J2=max(IM,JM)；
5. 将FLAG(J2)置为.FALSE.（j2 并入 i2，不再参与后续凝聚），按 2.2 的公式更新新簇与所有活跃簇的距离；
6. 对最近邻恰好指向被合并簇的簇重算最近邻；
7. NCL=NCL-1，若NCL>1则回到第 3 步，直到完成 N-1 次凝聚。

该实现同时输出完整的合并历史IA/IB/CRIT，这正是构建聚类树（第四节）所需的原始材料。

三、C++ 桥接层：Clustering 扩展模块（Clustering.cpp）

Clustering.cpp用 Boost.Python（nanobind 版本为nbClustering.cpp）导出两个函数，二者在构建配置 Code/ML/Cluster/Murtagh/CMakeLists.txt 中由rdkit_python_extension/rdkit_nanobind_extension编译为ML/Cluster目标下的同名扩展模块，并链接hc（由hc.c、hcdriver.c两个 f2c 转译文件构成）与RDGeneral库。

3.1 MurtaghCluster：从原始数据自动计算距离矩阵

static PyObject *Clustering_MurtaghCluster(python::object data, int nPts, int sz, int option);

• data：二维numpy.ndarray，形状为nPts × sz，即 nPts 个 sz 维数据点；
• 封装函数clusterit在 C++ 侧自行生成距离矩阵，而不是交给 f2c 转译的 Fortran 代码处理——源码注释直言这是为了避免"淹没在 f2c 的汪洋中"；
• 距离定义为欧氏距离平方：tmp = dataP[iTab+k] - dataP[jTab+k]; dists[pos] += tmp*tmp;，矩阵长度为n*(n-1)/2；
• 随后调用distdriver_(&n, &len, dists, &iopt, ia, ib, crit)进入 Fortran 核心。

3.2 MurtaghDistCluster：直接消费距离矩阵

static PyObject *Clustering_MurtaghDistCluster(python::object data, int nPts, int option);

• data：一维numpy.ndarray，长度必须为nPts*(nPts-1)/2，即按 2.1 节压缩格式存储的相异度向量；
• 该函数不再计算距离，直接调用distdriver_。这使上层可以在 Python 中自由选择任意距离度量（如基于分子指纹的 Tanimoto/Dice 相异度）后再进行层次聚类。

3.3 返回约定与内存管理

两个函数均返回一个三元素PyTuple：

• 元素 0/1：NPY_LONG数组ia、ib（合并的两个簇编号，1 基索引，长度 nPts，仅前 n-1 个有效）；
• 元素 2：NPY_DOUBLE数组crit（每次合并的相异度值）。

源码注释特别说明了内存策略：三个数组通过PyCapsule_New与PyArray_SetBaseObject挂接到 NumPy 数组的 base 对象上，由 Python 的引用计数统一释放，避免手动free与内存泄漏。

四、Python 接口层：Murtagh.py 与 ClusterData

rdkit/ML/Cluster/Murtagh.py 是对该扩展模块的高层封装，导入失败时（try/except ImportError）将MurtaghCluster/MurtaghDistCluster置为None，保证纯 Python 环境下的优雅降级——UnitTestCluster.py中也用unittest.skipIf对应处理。

4.1 方法常量与注册表

WARDS = 1 SLINK = 2 CLINK = 3 UPGMA = 4 MCQUITTY = 5 GOWER = 6 CENTROID = 7 methods = [ ("Ward's Minimum Variance", WARDS, "Ward's Minimum Variance"), ('Average Linkage', UPGMA, 'Group Average Linkage (UPGMA)'), ('Single Linkage', SLINK, 'Single Linkage (SLINK)'), ('Complete Linkage', CLINK, 'Complete Linkage (CLINK)'), # ("McQuitty",MCQUITTY,"McQuitty's method"), # ("Gower",GOWER,"Gower's median method"), ("Centroid", CENTROID, "Centroid method"), ]

注意：Fortran 层支持全部 7 种准则，但 Python 的methods注册表中 McQuitty 与 Gower 被注释禁用，官方对外推荐的方法是 Ward、UPGMA、SLINK、CLINK 与 Centroid 五种。

4.2 主入口 ClusterData

def ClusterData(data, nPts, method, isDistData=0):

参数说明：

函数内部先numpy.array(data)归一化，然后：

• 非距离数据：sz = data.shape[1]，调用MurtaghCluster(data, nPts, sz, method)；
• 距离数据：调用MurtaghDistCluster(data, nPts, method)；
• 最后用_ToClusters将ia/ib/crit转成聚类树。

4.3 结果转换：_ToClusters 与 _LookupDist

_LookupDist(dists, i, j, n)实现了 2.1 节的压缩索引读取：i==j返回 0，否则归一化i<j后返回dists[j*(j-1)//2 + i]。

_ToClusters则是合并历史的"树化"过程：

1. 为每个原始点创建叶子簇Clusters.Cluster(metric=0.0, data=i, index=i+1)；
2. 依次读取ia/ib/crit中记录的第 i 次合并（Fortran 数组为 1 基索引，Python 侧减 1），将cs[idx1]与cs[idx2]作为两个子节点，生成新簇Clusters.Cluster(metric=crit[i], children=[c1, c2], index=nPts+i+1)；
3. 返回单元素列表[newClust]，即整棵聚类树的根节点。

五、聚类树数据结构：Clusters.Cluster

rdkit/ML/Cluster/Clusters.py 定义了代表层次聚类树的Cluster类，是 Clustering 模块所有结果的承载容器。

5.1 属性一览

5.2 核心方法

• AddChild(child)/AddChildren(children)：添加子簇并触发_GenPoints()与_UpdateLength()；
• GetPoints()/GetPointsPositions()：惰性递归生成（首次调用后缓存于_points）；
• FindSubtree(index)：按索引深度优先查找子树；
• IsTerminal()：_len <= 1时为叶子节点；len(cluster)由__len__返回该簇覆盖的原始点数；
• GetChildren()：返回按Metric升序排序的子簇；
• Print(level, showData, offset)：以缩进树形打印簇结构；
• Compare(other, ignoreExtras)与__cmp__：提供比较能力，CMPTOL = 1e-6用于 metric 的浮点容差比较。

ClusterUtils.py中的GetNodeList返回按子树长度降序排列的全部节点列表，UnitTestCluster.py用一棵手工构造的树验证其顺序（[30, 31, 32, 10, 40, 41, 11, 1]）。

六、配套模块：距离矩阵、切分、可视化与标准化

Clustering 模块并不是孤立存在的，它依赖rdkit.ML.Cluster包中的若干配套模块（见 rdkit/ML/Cluster/）：

6.1 Resemblance.py：相异度矩阵的生成与索引约定

rdkit/ML/Cluster/Resemblance.py 文档字符串给出了压缩矩阵的完整约定：矩阵元素(row,col)（col>row）的线性索引为(col*(col-1))/2 + row，并给出了示例矩阵布局。提供的函数：

• EuclideanDistance(inData)：计算欧氏距离矩阵（注意返回的是欧氏距离即开根号后的值，而 C++ 侧MurtaghCluster计算的是距离平方，二者使用场景不同）；
• CalcMetricMatrix(inData, metricFunc)：对带有GetPosition()的对象列表生成度量矩阵；
• FindMinValInList(mat, nObjs)：在度量矩阵中查找最小值及其 (row, col)；
• ShowMetricMat(metricMat, nObjs)：以三角阵形式打印度量矩阵。

6.2 ClusterUtils.py：从聚类树切出 N 个簇

rdkit/ML/Cluster/ClusterUtils.py 除GetNodeList外，还提供：

• SplitIntoNClusters(cluster, n, breadthFirst=True)：把一棵聚类树切割为 n 个分支簇。breadthFirst=True用广度优先（_BreadthFirstSplit，每次选 metric 最大的簇按子节点 metric 降序裂开）；breadthFirst=False用高度优先（_HeightFirstSplit，每轮都重新按 metric 排序选最大者裂开）；
• FindClusterCentroidFromDists(cluster, dists)：找出簇内到其他点欧氏距离之和最小的点作为质心，并为每个点记录_distToCenter与_clustCenter。

UnitTestCluster.py的testDivide对一棵 5 簇树分别验证了两种切分方式的结果顺序，以及n超过簇大小的ValueError、n等于簇大小时返回全部叶子等边界行为。

6.3 ClusterVis.py：聚类树可视化

rdkit/ML/Cluster/ClusterVis.py 基于 Sping 绘图库将聚类树渲染为图。核心入口：

DrawClusterTree(cluster, canvas, size, ptColors=[], lineWidth=None, showIndices=0, showNodes=1, stopAtCentroids=0, logScale=0, tooClose=-1) ClusterToPDF(cluster, fileName, size=(300,300), ...) # 输出 PDF ClusterToSVG(cluster, fileName, size=(300,300), ...) # 输出 SVG ClusterToImg(cluster, fileName, size=(300,300), ...) # 输出位图（格式由文件名扩展名决定）

VisOpts类集中管理线宽、颜色、节点半径等绘制选项；logScale参数可对 metric 做对数缩放以改善树形观感（_scaleMetric实现log(metric^power / min)）。

6.4 Standardize.py：聚类前数据标准化

rdkit/ML/Cluster/Standardize.py 提供两种预处理方案：methods注册表包含"None"（不做标准化）与"Standard Deviation"（标准差标准化，委托rdkit.ML.Data.Stats.StandardizeMatrix实现），可在聚类前消除不同维度量纲差异的影响。

6.5 Butina.py：另一种聚类入口

rdkit/ML/Cluster/Butina.py 实现了 Butina（JCICS 39 747-750, 1999）的基于距离阈值的聚类算法，与 Murtagh 层次聚类互补：

def ClusterData(data, nPts, distThresh, isDistData=False, distFunc=EuclideanDist, reordering=False):

• distThresh：两点距离小于等于该值即视为邻居；
• 算法按邻居数降序挑选种子点，将种子及其所有未见邻居划为一个簇；
• reordering=True时，每生成一个簇就为剩余未分配分子重算未分配邻居数并重排，确保总是选择"未分配邻居最多"的分子作为下一个簇中心（UnitTestButina.py的test8_reordering_changes展示了开启重排序后簇划分从 3 簇变为 2 簇的行为差异）；
• 返回元组元组，每个簇的第一个元素是质心。

七、实战：基于分子指纹的层次聚类

将 Clustering 模块用于化合物集聚类时，仓库提供了现成的实用工具 rdkit/Chem/Fingerprints/ClusterMols.py。

7.1 核心调用链

from rdkit import DataStructs from rdkit.Chem.Fingerprints import FingerprintMols, MolSimilarity from rdkit.ML.Cluster import Murtagh dMat = GetDistanceMatrix(data, metric) # 指纹 → 压缩距离矩阵 clustTree = Murtagh.ClusterData(dMat, len(data), algorithmId, isDistData=1)[0]

• GetDistanceMatrix将指纹两两比较得到相异度矩阵：相似度度量下取1.0 - metric(fp1, fp2)；若两个指纹位数不同，会先调用DataStructs.FoldFingerprint折叠到相同位数再比较；
• 距离矩阵按Resemblance.py的压缩格式存储，随后交给Murtagh.ClusterData(..., isDistData=1)进行层次聚类；
• 结果树上每个叶节点通过SetName打上分子标签，有活性数据时还会SetData存储活性值，供可视化时着色。

7.2 命令行用法

ClusterMols.py支持从数据库或文本文件读入指纹数据直接聚类，常用参数：

典型示例（文件头部注释给出）：

python ClusterMols.py -d data.gdb -t daylight_sig --idName="CAS_TF" \ -o clust1.pkl --actTable="dop_test" --actName="moa_quant"

7.3 结果树的使用

拿到聚类树根节点后，常用组合操作：

from rdkit.ML.Cluster import ClusterUtils # 打印整棵树 clustTree.Print(showData=True) # 切成 N 个簇 clusters = ClusterUtils.SplitIntoNClusters(clustTree, 10) # 计算各簇质心（需提供距离矩阵） centroid = ClusterUtils.FindClusterCentroidFromDists(clustTree, dMat)

八、测试验证：行为契约一览

仓库用多组测试锁定 Clustering 模块的行为：

• rdkit/ML/Cluster/UnitTestCluster.py：使用 Romesburg《Cluster Analysis for Researchers》一书中的经典 5 点数据集（[[10,5],[20,20],[30,10],[30,15],[5,10]]）验证Murtagh.ClusterData三种输入方式（原始数据、距离矩阵、以及同一数据不同调用路径）得到完全等价的聚类树（通过Compare(ignoreExtras=0)断言）；
• rdkit/ML/Cluster/murtagh_test.py：对同一数据集分别以WARDS、SLINK、CLINK、UPGMA运行并Print()输出，可直观对比不同准则下的树形；
• rdkit/ML/Cluster/UnitTestButina.py：覆盖 Butina 聚类的多个边界场景，包括全单例、全一簇、零距离、空输入，以及数据维度不匹配时的ValueError、非法类型时的TypeError；
• 测试数据目录 rdkit/ML/Cluster/test_data/ 中存放了ferro.UPGMA.pkl、struct.Wards.pkl等预期聚类树 pickle，用于回归比对。

九、构建与依赖说明

扩展模块的构建定义位于 Code/ML/Cluster/Murtagh/CMakeLists.txt：

• rdkit_library(hc hc.c hcdriver.c)：将 f2c 转译的 Fortran 代码编译为静态库hc（hc.f经 f2c 生成hc.c，hcdriver.f生成hcdriver.c）；
• 当开启RDK_BUILD_BOOST_PYTHON_WRAPPERS时编译 Boost.Python 版Clustering扩展；开启RDK_BUILD_NANOBIND_WRAPPERS时编译 nanobind 版nbClustering.cpp；
• 扩展链接hc与RDGeneral，安装目标为ML/Cluster目录，即 Python 侧的rdkit.ML.Cluster.Clustering。

这也解释了Murtagh.py中try/except的必要性：如果构建时未启用 Python 包装器，或hc库未编译，则MurtaghCluster/MurtaghDistCluster为None，ClusterData将不可用，调用方需自行降级处理。

十、小结

rdkit.ML.Cluster.Clustering是 RDKit 层次聚类能力的底层引擎：Fortran 的HC凝聚算法提供 Ward/SLINK/CLINK/UPGMA 等七种准则与完整的合并历史，C++ 桥接层负责距离矩阵生成与 NumPy 数据交换，Python 层的Murtagh.ClusterData则把合并历史转成易于操作的Cluster聚类树。理解这条链路后，读者不仅可以正确选用isDistData、压缩距离矩阵等关键约定，还能借助ClusterUtils.SplitIntoNClusters、ClusterVis与Butina等配套模块，在分子指纹相似性聚类场景中构建完整可用的工作流。

• 科学计算
• 科研
• 机器学习