基于SMARTS从SMILES检测化学基团的技术问题求助

现状与需求 #

现有基于RDKit的代码可有效检测常规化学基团，但存在以下局限性：

• 无法识别芳香环、杂环、脂环（如环己醇环）等环结构
• 无法区分烯烃的顺反构型与芳香烃中的双键
• 需要实现：按环大小、是否含杂原子、是否芳香分类检测环结构；准确判定烯烃的顺反构型

现有代码 #

from rdkit import Chem

def find_smiles_patterns(smiles):
    mol = Chem.MolFromSmiles(smiles)
    if mol is None:
        return "Invalid SMILES string. Unable to parse molecule."

    # 存储检测到的化学基团
    chemical_groups = []

    # 用于识别基团的SMARTS模式
    smarts_patterns = {
        'C=C': 'Alkene',
        '[CX2]#[CX2]': 'Alkyne',
        '[CX3]=[CX2]=[CX3]': 'Allene',
        '[ClX1][CX4]': 'Alkylchloride',
        '[FX1][CX4]': 'Alkylfluoride',
        '[BrX1][CX4]': 'Alkylbromide',
        '[IX1][CX4]': 'Alkyliodide',
        '[OX2H][CX4H2;!$(C([OX2H])[O,S,#7,#15])]': 'Primary_alcohol',
        '[OX2H][CX4H;!$(C([OX2H])[O,S,#7,#15])]': 'Secondary_alcohol',
        '[OX2H][CX4D4;!$(C([OX2H])[O,S,#7,#15])]': 'Tertiary_alcohol',
        '[OX2]([CX4;!$(C([OX2])[O,S,#7,#15,F,Cl,Br,I])])[CX4;!$(C([OX2])[O,S,#7,#15])]': 'Dialkylether',
        '[SX2]([CX4;!$(C([OX2])[O,S,#7,#15,F,Cl,Br,I])])[CX4;!$(C([OX2])[O,S,#7,#15])]': 'Dialkylthioether',
        '[OX2](c)[CX4;!$(C([OX2])[O,S,#7,#15,F,Cl,Br,I])]': 'Alkylarylether',
        '[c][OX2][c]': 'Diarylether',
        '[SX2](c)[CX4;!$(C([OX2])[O,S,#7,#15,F,Cl,Br,I])]': 'Alkylarylthioether',
        '[c][SX2][c]': 'Diarylthioether',
        '[O+;!$([O]~[!#6]);!$([S]*~[#7,#8,#15,#16])]': 'Oxonium',
        '[NX3H2+0,NX4H3+;!$([N][!C]);!$([N]*~[#7,#8,#15,#16])]': 'Primary_aliph_amine',
        '[NX3H1+0,NX4H2+;!$([N][!C]);!$([N]*~[#7,#8,#15,#16])]': 'Secondary_aliph_amine',
        '[NX3H0+0,NX4H1+;!$([N][!C]);!$([N]*~[#7,#8,#15,#16])]': 'Tertiary_aliph_amine',
        '[NX4H0+;!$([N][!C]);!$([N]*~[#7,#8,#15,#16])]': 'Quaternary_aliph_ammonium',
        '[!#6;!R0]': 'Heterocyclic'
        # 可扩展更多模式
    }

    # 基于IUPAC命名的基团检测优先级
    priority_order = [
        'Carboxylic_acid',
        'Carboxylic_ester',
        'Lactone',
        'Carboxylic_anhydride',
        'Carbothioic_acid',
        'Aldehyde',
        'Ketone',
        'Alkylchloride',
        'Alkylfluoride',
        'Alkylbromide',
        'Alkyliodide',
        'Alcohol',
        'Primary_alcohol',
        'Secondary_alcohol',
        'Tertiary_alcohol',
        'Dialkylether',
        'Alkene',
        'Alkyne',
        'Allene',
        'Dialkylthioether',
        'Alkylarylether',
        'Diarylether',
        'Alkylarylthioether',
        'Diarylthioether',
        'Oxonium',
        'Primary_aliph_amine',
        'Secondary_aliph_amine',
        'Tertiary_aliph_amine',
        'Quaternary_aliph_ammonium',
        'Heterocycle'
        # 可扩展更多优先级
    ]

    # 跟踪已检测的原子索引，避免重复计数
    atom_indices = set()
    
    # 按优先级遍历检测基团
    for group in priority_order:
        if group in smarts_patterns.values():
            for smarts_pattern, chemical_group in smarts_patterns.items():
                if chemical_group == group:
                    pattern = Chem.MolFromSmarts(smarts_pattern)
                    if pattern:
                        matches = mol.GetSubstructMatches(pattern)
                        if len(matches) > 0:
                            print('matches !!! : ', smarts_pattern , smarts_patterns[smarts_pattern], pattern)
                            print(matches,'\n\n')
                            for match in matches:
                                match_set = set(match)
                                if not any(atom_index in match_set for atom_index in atom_indices):
                                    chemical_groups.append(chemical_group)
                                    atom_indices.update(match_set)

    return chemical_groups

smiles = "c1(cccc2c1ccc1c2cccc1CCl)CO"
print(find_smiles_patterns(smiles))

测试结果 #

matches !!! :  [ClX1][CX4] Alkylchloride <rdkit.Chem.rdchem.Mol object at .............>
((15, 14),) 

matches !!! :  [OX2H][CX4H2;!$(C([OX2H])[O,S,#7,#15])] Primary_alcohol <rdkit.Chem.rdchem.Mol object at .............>
((17, 16),) 

['Alkylchloride', 'Primary_alcohol']

解决方案 #

1. 环结构检测的SMARTS模式 #

以下是按分类整理的SMARTS，可直接添加到smarts_patterns字典中：

芳香环类

• 六元芳香碳环（苯环）：c1ccccc1 → 对应名称Benzene_ring
• 五元芳香碳环：c1cccc1 → 对应名称Cyclopentadienyl_ring
• 六元含氮杂芳香环（吡啶环）：c1ccncc1 → 对应名称Pyridine_ring
• 五元含氮杂芳香环（吡咯环）：c1ccnc1 → 对应名称Pyrrole_ring
• 通用芳香环（任意大小，含杂原子）：[r]1~[r]~[r]~[r]~[r]1（五元）、[r]1~[r]~[r]~[r]~[r]~[r]1（六元）→ 对应名称Aromatic_ring

脂环类

• 六元脂环碳环（环己烷环）：C1CCCCC1 → 对应名称Cyclohexane_ring
• 五元脂环碳环（环戊烷环）：C1CCCC1 → 对应名称Cyclopentane_ring
• 六元含氮杂脂环（哌啶环）：C1CCNCC1 → 对应名称Piperidine_ring
• 通用脂环（任意大小，含杂原子）：[R;!r]1~[R;!r]~[R;!r]~[R;!r]~[R;!r]1（五元）、[R;!r]1~[R;!r]~[R;!r]~[R;!r]~[R;!r]~[R;!r]1（六元）→ 对应名称Aliphatic_cyclic_ring

按环大小区分

可通过[R<n>]指定环原子数，比如[R6]表示六元环原子，因此六元芳香碳环也可写为：[cR6]1[cR6][cR6][cR6][cR6][cR6]1

2. 烯烃顺反构型判定方案 #

RDKit可通过解析双键的立体化学信息区分顺反，需在代码中添加以下逻辑：

def detect_alkene_stereochemistry(mol):
    # 先分配立体化学信息
    Chem.AssignStereochemistry(mol, cleanIt=True, force=True)
    stereochem_groups = []
    for bond in mol.GetBonds():
        # 只处理非芳香的双键
        if bond.GetBondType() == Chem.BondType.DOUBLE and not bond.GetIsAromatic():
            stereo = bond.GetStereo()
            if stereo == Chem.BondStereo.STEREOZ:
                stereochem_groups.append('Alkene_Z')
            elif stereo == Chem.BondStereo.STEREOE:
                stereochem_groups.append('Alkene_E')
    return stereochem_groups

然后在find_smiles_patterns函数中调用该方法，将结果合并到chemical_groups中（注意需移除原有普通烯烃的检测逻辑，避免重复）。

代码整合注意事项 #

• 将环结构的SMARTS添加到smarts_patterns后，需同步更新priority_order列表，调整检测优先级（比如环结构优先级低于官能团）
• 环结构检测时，需注意避免重复计数同一个环，可通过跟踪环的ID（mol.GetRingInfo()）来实现，而非仅跟踪原子索引

内容的提问来源于stack exchange，提问作者Cedric Baerlocher