去盐与去重复：化合物数据清洗的第一步

拿到一批化合物数据，第一件事是去盐与去重。这两步看似简单，但「怎么去重」「重复的活性数据怎么合并」这些判断，直接影响后续模型的质量。

去盐：保留哪个片段

from rdkit import Chem
from rdkit.Chem.MolStandardize import rdMolStandardize

# ---- 方法一：保留最大的有机片段（推荐）----
def desalt(smiles):
    mol = Chem.MolFromSmiles(smiles)
    if mol is None:
        return None
    parent = rdMolStandardize.FragmentParent(mol)
    return Chem.MolToSmiles(parent) if parent else None

print(desalt("CC(=O)[O-].[Na+]"))          # CC(=O)[O-]
print(desalt("CN1CCC[C@H]1c1cccnc1.Cl"))   # 去掉 HCl

# ---- 方法二：用盐列表移除 ----
remover = rdMolStandardize.LargestFragmentChooser()
mol = remover.choose(Chem.MolFromSmiles("CC(=O)O.O.O"))

# ---- 需要注意的情况 ----
#
# 【情况一：共晶与共药】
#   两个都是活性成分的共晶
#   → 【去掉一个就丢失了信息】
#   → 应该标记出来单独处理
#
# 【情况二：季铵盐等本身带电的化合物】
#   CC[N+](C)(C)C.[Br-]
#   → 阳离子部分是「主体」，溴离子是抗衡离子
#   → FragmentParent 会正确处理
#
# 【情况三：两个片段大小相近】
#   「最大片段」的判断可能不符合预期
#   → 【应该检查】

def check_desalting(smiles_list, n_sample=20):
    """检查去盐是否合理"""
    multi_component = []
    for smi in smiles_list:
        if "." not in smi:
            continue
        mol = Chem.MolFromSmiles(smi)
        if mol is None:
            continue
        frags = Chem.GetMolFrags(mol, asMols=True)
        sizes = sorted([f.GetNumHeavyAtoms() for f in frags], reverse=True)
        # 【如果前两大片段大小接近，需要人工检查】
        if len(sizes) >= 2 and sizes[1] > 0.5 * sizes[0]:
            multi_component.append((smi, sizes))
    print(f"【{len(multi_component)} 个分子的多个片段大小接近，需检查】")
    return multi_component[:n_sample]

去重：用什么作为标识

import pandas as pd
from rdkit import Chem
from rdkit.Chem.MolStandardize import rdMolStandardize

def make_key(smiles, level="standard"):
    """生成不同严格程度的去重键"""
    mol = Chem.MolFromSmiles(smiles)
    if mol is None:
        return None
    mol = rdMolStandardize.Cleanup(mol)
    mol = rdMolStandardize.FragmentParent(mol)
    if mol is None:
        return None

    if level == "standard":
        return Chem.MolToSmiles(mol)
    if level == "no_stereo":
        Chem.RemoveStereochemistry(mol)
        return Chem.MolToSmiles(mol)
    if level == "inchikey":
        return Chem.InchiToInchiKey(Chem.MolToInchi(mol))
    if level == "skeleton":
        return Chem.InchiToInchiKey(Chem.MolToInchi(mol))[:14]
    raise ValueError(level)

# 【比较不同严格程度下的重复数】
def compare_dedup_levels(smiles_list):
    for level in ["standard", "no_stereo", "inchikey", "skeleton"]:
        keys = [make_key(s, level) for s in smiles_list]
        keys = [k for k in keys if k]
        print(f"{level:12s}: {len(keys)} → {len(set(keys))} "
              f"（去掉 {len(keys) - len(set(keys))} 个重复）")

# 【判断该用哪个级别】：
#   如果数据中立体化学标注不完整
#   → 用 no_stereo 或 skeleton，避免同一化合物
#     因为标注差异被分成两条
#   如果立体化学与活性相关且标注完整
#   → 用 standard 或 inchikey

重复数据的活性合并：关键判断

# 【发现重复后，活性值怎么处理？】
#
# 【第一步：先检查一致性，不要直接合并】

import numpy as np

def check_duplicate_consistency(df, key_col="std_smiles",
                                value_col="pActivity", threshold=1.0):
    """检查重复测定的一致性"""
    grouped = df.groupby(key_col)[value_col]
    stats = grouped.agg(["count", "min", "max", "std", "median"])
    stats["spread"] = stats["max"] - stats["min"]

    dups = stats[stats["count"] > 1]
    print(f"有重复的分子: {len(dups)}")
    print(f"重复次数分布:\n{dups['count'].value_counts().head()}")

    inconsistent = dups[dups["spread"] > threshold]
    print(f"\n【差异 > {threshold} 个 log 单位的: {len(inconsistent)}】")
    print(inconsistent.sort_values("spread", ascending=False).head(10))

    return stats, inconsistent

# 【判读】：
#   spread < 0.3   → 一致性好，可以合并
#   spread 0.3~1.0 → 正常的实验波动
#   spread > 1.0   → 【可能有问题】：
#     - 不同的测定方法/条件
#     - 数据录入错误
#     - 【结构标识错误】（两个不同的化合物被当成同一个）
#     - 化合物不稳定
#   → 【应该单独检查，考虑剔除或分开处理】

# 【第二步：选择合并策略】
#
# 中位数（推荐）：
#   对异常值稳健
#   df.groupby(key)[value].median()
#
# 均值：
#   受异常值影响
#
# 最大值：
#   假设「测不到活性是实验问题」
#   → 【会系统性高估】
#
# 保留信息最多的那条：
#   如按数据来源的可靠性排序
#
# 【最好的做法】：
#   合并后【保留重复次数与离散程度】作为额外列
#   → 建模时可以给一致性好的数据更高权重
#   → 【也可以用离散程度估计实验噪声的下限】

def merge_duplicates(df, key_col="std_smiles", value_col="pActivity"):
    agg = df.groupby(key_col).agg(
        value=(value_col, "median"),
        n_measurements=(value_col, "count"),
        spread=(value_col, lambda x: x.max() - x.min()),
    ).reset_index()
    return agg

# 【实验噪声的估计】：
#   重复测定的标准差 = 实验噪声的下限
#   → 【模型的 RMSE 不可能低于这个值】
#   → 这给了你一个合理的性能上限参照（见 044）

跨数据源合并时的额外问题

• 不同来源的实验条件可能不同：即使化合物相同，IC50 也可能差几倍（见 156《AI ADMET 预测模型》）。合并前应该检查各来源的分布是否一致；
• 标识符的映射：CAS 号可能对应盐型，数据库 ID 之间的映射不完全（见 033《InChI 与 InChIKey》）；
• 数据来源作为特征：如果无法统一，把「来源」作为一个特征，让模型显式建模体系差异；
• 保留追溯信息：合并后记录每条数据来自哪些原始记录——出问题时能查；
• 不要过早合并：先分别分析各来源，确认可比后再合并。

完整的清洗流程

def clean_dataset(df, smiles_col="smiles", value_col="pActivity"):
    """完整的数据清洗流程"""
    n0 = len(df)
    log = {}

    # 1) 解析与标准化
    df["std_smiles"] = df[smiles_col].apply(lambda s: make_key(s, "standard"))
    df = df.dropna(subset=["std_smiles"])
    log["解析/标准化失败"] = n0 - len(df)

    # 2) 【检查异常结构】（见 048）
    #    分子量范围、原子类型、片段数等
    from rdkit.Chem import Descriptors
    def is_reasonable(smi):
        m = Chem.MolFromSmiles(smi)
        if m is None:
            return False
        mw = Descriptors.MolWt(m)
        return 50 < mw < 2000 and m.GetNumHeavyAtoms() >= 3
    n1 = len(df)
    df = df[df["std_smiles"].apply(is_reasonable)]
    log["异常结构"] = n1 - len(df)

    # 3) 【检查重复一致性】
    stats, inconsistent = check_duplicate_consistency(df)
    #    可选：剔除高度不一致的
    bad_keys = set(inconsistent.index)
    n2 = len(df)
    df = df[~df["std_smiles"].isin(bad_keys)]
    log["不一致的重复"] = n2 - len(df)

    # 4) 合并重复
    df = merge_duplicates(df, value_col=value_col)
    log["合并后唯一分子数"] = len(df)

    # 5) 【记录清洗日志】
    print("清洗日志:")
    for k, v in log.items():
        print(f"  {k}: {v}")
    return df

# 【清洗日志的价值】：
#   1) 可复现
#   2) 【如果丢失的数据比例很高，说明数据源有问题】
#   3) 出问题时能定位

延伸资源

• 分子标准化：046《分子标准化》；异常结构：048《异常结构识别》；InChI：033《InChI 与 InChIKey》；
• QSAR：044《QSAR 入门》；ADMET 数据：156《AI ADMET 预测模型》。