用 Tanimoto 相似性做化合物聚类

化合物聚类有三个典型用途：理解库的结构组成、去冗余、挑选多样化子集。Butina 聚类是化学信息学中最常用的方法——它不需要预设簇数，且对化学相似性的定义直观。

完整流程

from rdkit import Chem, DataStructs
from rdkit.Chem import rdFingerprintGenerator
from rdkit.ML.Cluster import Butina
import numpy as np

def butina_cluster(smiles_list, cutoff=0.35, radius=2, fp_size=2048):
    """cutoff 是距离阈值：距离 = 1 - Tanimoto 相似度
       cutoff=0.35 ≈ 相似度 0.65 以上算同簇"""
    gen = rdFingerprintGenerator.GetMorganGenerator(radius=radius, fpSize=fp_size)
    mols = [Chem.MolFromSmiles(s) for s in smiles_list]
    valid_idx = [i for i, m in enumerate(mols) if m is not None]
    fps = [gen.GetFingerprint(mols[i]) for i in valid_idx]

    # 构建下三角距离矩阵（Butina 要求的格式）
    n = len(fps)
    dists = []
    for i in range(1, n):
        sims = DataStructs.BulkTanimotoSimilarity(fps[i], fps[:i])
        dists.extend(1.0 - s for s in sims)

    clusters = Butina.ClusterData(dists, n, cutoff, isDistData=True)
    # 返回：每个簇是一个索引元组，第一个元素是簇中心
    return clusters, valid_idx

clusters, valid_idx = butina_cluster(df["clean_smiles"].tolist(), cutoff=0.35)
print(f"{len(valid_idx)} 个分子 → {len(clusters)} 个簇")
print(f"最大簇 {len(clusters[0])} 个成员；单例簇 "
      f"{sum(1 for c in clusters if len(c)==1)} 个")

Butina 算法在做什么

• 对每个分子，统计与它距离小于 cutoff 的邻居数；
• 取邻居最多的分子作为第一个簇中心，把它和它的所有邻居归为一簇；
• 从剩余分子中重复，直到全部分配完；
• 特点：确定性（无随机性）、不需预设簇数、每个簇有明确的中心分子。

阈值怎么选

cutoff（距离）	对应相似度	效果	适用
0.2	> 0.8	簇很多很小，只有近似类似物同簇	精细去重
0.35	> 0.65	常用平衡点	多数场景
0.5	> 0.5	簇较大较粗	粗略分组
0.7	> 0.3	大部分分子挤在少数簇	基本无用

阈值强烈依赖指纹类型——同样的 0.35，用 ECFP4 和用 MACCS 键得到的簇结构完全不同。换指纹就要重新校准阈值。实践建议：试 3～4 个值，看簇数分布是否合理（既不是一个巨簇，也不是全是单例）。

三个实际用途

# 用途一：挑代表分子（去冗余）
representatives = [c[0] for c in clusters]     # 每簇第一个是中心
print(f"从 {len(valid_idx)} 个降到 {len(representatives)} 个代表")
rep_smiles = [df["clean_smiles"].iloc[valid_idx[i]] for i in representatives]

# 用途二：多样化采样（预算有限时买哪些）
# 按簇大小排序，优先从大簇取（覆盖主要化学型），
# 再补单例簇（覆盖独特结构）
sorted_clusters = sorted(clusters, key=len, reverse=True)
budget = 30
picked = []
for c in sorted_clusters:
    if len(picked) >= budget:
        break
    picked.append(c[0])
print(f"多样化挑选 {len(picked)} 个，覆盖 {len(picked)} 个不同化学型")

# 用途三：理解库的组成
sizes = [len(c) for c in clusters]
print(f"簇大小分布：中位 {np.median(sizes):.0f}，最大 {max(sizes)}，"
      f"单例占比 {sum(1 for s in sizes if s==1)/len(sizes):.1%}")
# 单例占比很高 → 库很多样
# 少数巨簇 → 库中有大量类似物，可能来自同一系列

大库的处理策略

Butina 需要完整的两两距离矩阵，复杂度是 O(n²)。超过约 5 万个分子就会内存吃紧：

# 策略一：分批聚类再合并
# 先按某个粗特征（如骨架、分子量段）分组，组内聚类

# 策略二：用 MiniBatchKMeans 等可扩展方法
from sklearn.cluster import MiniBatchKMeans
X = np.array([np.array(fp) for fp in fps])         # 转成数组
km = MiniBatchKMeans(n_clusters=1000, batch_size=1000, random_state=0)
labels = km.fit_predict(X)
# 注意：KMeans 用欧氏距离，与 Tanimoto 的化学含义不同，
#       但对大规模粗分组够用

# 策略三：先用 MaxMin 挑多样化子集，再聚类
from rdkit.SimDivFilters import MaxMinPicker
picker = MaxMinPicker()
idx = picker.LazyBitVectorPick(fps, len(fps), 5000, seed=42)
# MaxMin 直接给出最多样化的 N 个分子，不需要完整距离矩阵

# 策略四：用骨架分组（最快的粗分组）
from rdkit.Chem.Scaffolds import MurckoScaffold
df["scaffold"] = df["clean_smiles"].map(
    lambda s: MurckoScaffold.MurckoScaffoldSmiles(s) or "")
print(f"{df['scaffold'].nunique()} 个不同骨架")

MaxMinPicker 是被低估的工具：直接挑出最多样化的 N 个分子，不需要完整距离矩阵，非常适合「预算有限、要买最多样的 30 个」这类需求。

注意事项

• 聚类是探索手段，不是硬分类：簇边界是阈值决定的，不代表本质区别。
• 结构相似 ≠ 活性相似：同一簇内可能有活性悬崖。
• 先标准化再聚类：否则盐型不同的同一化合物会被分到不同簇。
• 用于数据划分时要谨慎：按簇划分能减少泄漏，但骨架划分（见 332《用 Scaffold Split 评估模型真实外推能力》）更常用也更标准。

延伸资源

• 概念：039《Tanimoto 相似性》、050《分子相似性聚类》；指纹：330《计算 ECFP4 Fingerprint》、036《Morgan Fingerprint 入门》；
• 数据划分：332《用 Scaffold Split 评估模型真实外推能力》。