05 - 复杂分子取代基替换¶

本 Notebook 使用真实的 Gaussian 16 过渡态输出 2-TS1-Sp.log 作为母体，演示悬挂取代基替换 API。每个结果都保留三维 conformer，并由 MolsToGridDofImage 直接从带构象的 rdmol 绘制为 SVG。

稳定的单点替换要求 query 描述一个与母体只有一条外部非环键相连的连通悬挂片段；最后还演示实验性的双锚点整体替换。下面覆盖 query 单点替换、bind_idx、replace_all=True、显式 start_idx / end_idx、金属配位键、双键 E / Z 立体策略和错误边界。

In [1]:

from pathlib import Path

import pandas as pd
from rdkit.Chem import AllChem
from rdkit.Chem.Draw import IPythonConsole
from IPython.display import SVG, display
from rdkit import Chem
from rdkit_dof import MolsToGridDofImage
import py3Dmol

from molop import AutoParser, molopconfig
from molop.structure import replace_substituent
from molop.structure.substitution_sites import find_sites, normalize_query

# 让 SVG 图中显示 RDKit 的 0-based 原子索引，便于对应 bind_idx/start_idx/end_idx。
IPythonConsole.drawOptions.addAtomIndices = True
molopconfig.quiet()
sample_candidates = (
    Path("tests/test_files/g16log/2-TS1-Sp.log"),
    Path("../../../tests/test_files/g16log/2-TS1-Sp.log"),
    Path("2-TS1-Sp.log"),
)
sample_path = next((path for path in sample_candidates if path.is_file()), None)
if sample_path is None:
    raise FileNotFoundError("请将 2-TS1-Sp.log 放在 Notebook 同目录，或从 MolOP 源码树运行。")

batch = AutoParser(sample_path, parser_detection="g16log", n_jobs=1, only_last_frame=True)
frame = batch[0][-1]
parent = frame.rdmol
if parent is None:
    raise ValueError("Gaussian frame 没有得到可编辑的 RDKit 分子。")
parent = Chem.Mol(parent)

def show_svg(molecules, legends):
    if any(
        molecule.GetNumConformers() == 0 or not molecule.GetConformer().Is3D()
        for molecule in molecules
    ):
        raise ValueError("SVG 演示要求所有分子都保留三维 conformer。")
    svg = MolsToGridDofImage(
        molecules,
        molsPerRow=min(2, len(molecules)),
        subImgSize=(720, 520),
        legends=legends,
        use_svg=True,
        return_image=False,
    )
    if not isinstance(svg, str):
        raise RuntimeError("rdkit-dof 没有返回 SVG 图像。")
    display(SVG(data=svg))

def show_3d(molecules, legends):
    if len(molecules) != len(legends):
        raise ValueError("molecules and legends must have the same length.")
    viewer = py3Dmol.view(
        width=720 * len(molecules),
        height=520,
        viewergrid=(1, len(molecules)),
        linked=False,
    )
    for column, (molecule, _legend) in enumerate(zip(molecules, legends, strict=True)):
        viewer.addModel(
            Chem.MolToMolBlock(molecule, confId=0),
            "mol",
            viewer=(0, column),
        )
        viewer.setStyle({"stick": {"radius": 0.16}}, viewer=(0, column))
        conformer = molecule.GetConformer()
        for atom in molecule.GetAtoms():
            position = conformer.GetAtomPosition(atom.GetIdx())
            viewer.addLabel(
                str(atom.GetIdx()),
                {
                    "position": {"x": position.x, "y": position.y, "z": position.z},
                    "fontSize": 8,
                    "fontColor": "black",
                    "backgroundColor": "white",
                    "backgroundOpacity": 0.65,
                    "showBackground": True,
                    "inFront": True,
                },
                viewer=(0, column),
            )
        viewer.zoomTo(viewer=(0, column))
    viewer.show()

print("输入文件：", sample_path)
print("母体 frame：", frame.frame_id)
print("母体原子数 / 键数 / 分子式：", parent.GetNumAtoms(), parent.GetNumBonds(), frame.formula)
print("母体 conformer 数 / 是否为 3D：", parent.GetNumConformers(), parent.GetConformer().Is3D())
show_svg([parent], ["2-TS1-Sp.log | parent"])
show_3d([parent], ["2-TS1-Sp.log | parent"])

from pathlib import Path import pandas as pd from rdkit.Chem import AllChem from rdkit.Chem.Draw import IPythonConsole from IPython.display import SVG, display from rdkit import Chem from rdkit_dof import MolsToGridDofImage import py3Dmol from molop import AutoParser, molopconfig from molop.structure import replace_substituent from molop.structure.substitution_sites import find_sites, normalize_query # 让 SVG 图中显示 RDKit 的 0-based 原子索引，便于对应 bind_idx/start_idx/end_idx。 IPythonConsole.drawOptions.addAtomIndices = True molopconfig.quiet() sample_candidates = ( Path("tests/test_files/g16log/2-TS1-Sp.log"), Path("../../../tests/test_files/g16log/2-TS1-Sp.log"), Path("2-TS1-Sp.log"), ) sample_path = next((path for path in sample_candidates if path.is_file()), None) if sample_path is None: raise FileNotFoundError("请将 2-TS1-Sp.log 放在 Notebook 同目录，或从 MolOP 源码树运行。") batch = AutoParser(sample_path, parser_detection="g16log", n_jobs=1, only_last_frame=True) frame = batch[0][-1] parent = frame.rdmol if parent is None: raise ValueError("Gaussian frame 没有得到可编辑的 RDKit 分子。") parent = Chem.Mol(parent) def show_svg(molecules, legends): if any( molecule.GetNumConformers() == 0 or not molecule.GetConformer().Is3D() for molecule in molecules ): raise ValueError("SVG 演示要求所有分子都保留三维 conformer。") svg = MolsToGridDofImage( molecules, molsPerRow=min(2, len(molecules)), subImgSize=(720, 520), legends=legends, use_svg=True, return_image=False, ) if not isinstance(svg, str): raise RuntimeError("rdkit-dof 没有返回 SVG 图像。") display(SVG(data=svg)) def show_3d(molecules, legends): if len(molecules) != len(legends): raise ValueError("molecules and legends must have the same length.") viewer = py3Dmol.view( width=720 * len(molecules), height=520, viewergrid=(1, len(molecules)), linked=False, ) for column, (molecule, _legend) in enumerate(zip(molecules, legends, strict=True)): viewer.addModel( Chem.MolToMolBlock(molecule, confId=0), "mol", viewer=(0, column), ) viewer.setStyle({"stick": {"radius": 0.16}}, viewer=(0, column)) conformer = molecule.GetConformer() for atom in molecule.GetAtoms(): position = conformer.GetAtomPosition(atom.GetIdx()) viewer.addLabel( str(atom.GetIdx()), { "position": {"x": position.x, "y": position.y, "z": position.z}, "fontSize": 8, "fontColor": "black", "backgroundColor": "white", "backgroundOpacity": 0.65, "showBackground": True, "inFront": True, }, viewer=(0, column), ) viewer.zoomTo(viewer=(0, column)) viewer.show() print("输入文件：", sample_path) print("母体 frame：", frame.frame_id) print("母体原子数 / 键数 / 分子式：", parent.GetNumAtoms(), parent.GetNumBonds(), frame.formula) print("母体 conformer 数 / 是否为 3D：", parent.GetNumConformers(), parent.GetConformer().Is3D()) show_svg([parent], ["2-TS1-Sp.log | parent"]) show_3d([parent], ["2-TS1-Sp.log | parent"])

输入文件： ../../../tests/test_files/g16log/2-TS1-Sp.log
母体 frame： 0
母体原子数 / 键数 / 分子式： 101 107 C44H45F3NNiO4PS2
母体 conformer 数 / 是否为 3D： 1 True

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

1. 发现真实悬挂取代基位点¶

母体中可以识别出磺酰基上的 CF₃、叔丁基，以及连接到膦或螯合碳的苯基。稠合芳基因为存在多条外部连接，不会被误判为可替换的悬挂片段。

In [2]:

def site_rows(query):
    query_mol = normalize_query(query)
    rows = []
    for site in find_sites(parent, query_mol):
        scaffold = parent.GetAtomWithIdx(site.scaffold_atom)
        substituent = parent.GetAtomWithIdx(site.substituent_atom)
        rows.append({
            "query": query,
            "scaffold_idx": site.scaffold_atom,
            "scaffold_element": scaffold.GetSymbol(),
            "substituent_idx": site.substituent_atom,
            "substituent_element": substituent.GetSymbol(),
            "bond": str(site.bond_type),
            "matched_atoms": len(site.query_match),
        })
    return rows

site_table = pd.DataFrame(
    site_rows("C(F)(F)F") + site_rows("C(C)(C)C") + site_rows("c1ccccc1")
)
display(site_table)

def site_rows(query): query_mol = normalize_query(query) rows = [] for site in find_sites(parent, query_mol): scaffold = parent.GetAtomWithIdx(site.scaffold_atom) substituent = parent.GetAtomWithIdx(site.substituent_atom) rows.append({ "query": query, "scaffold_idx": site.scaffold_atom, "scaffold_element": scaffold.GetSymbol(), "substituent_idx": site.substituent_atom, "substituent_element": substituent.GetSymbol(), "bond": str(site.bond_type), "matched_atoms": len(site.query_match), }) return rows site_table = pd.DataFrame( site_rows("C(F)(F)F") + site_rows("C(C)(C)C") + site_rows("c1ccccc1") ) display(site_table)

	query	scaffold_idx	scaffold_element	substituent_idx	substituent_element	bond	matched_atoms
0	C(F)(F)F	94	S	97	C	SINGLE	4
1	C(C)(C)C	32	S	44	C	SINGLE	4
2	c1ccccc1	3	P	10	C	SINGLE	6
3	c1ccccc1	3	P	21	C	SINGLE	6
4	c1ccccc1	5	C	33	C	SINGLE	6

2. Query 唯一匹配：CF₃ → 氧化取代基¶

replacement 使用 RDKit 支持的哑原子连接标记。[*]COC(=O)C 中的 [*] 表示与母体相连的位置；装配时会移除该标记而保留真实取代基。

In [3]:

cf3_replaced = replace_substituent(
    parent, query="C(F)(F)F", replacement="[*]COC(=O)C",
    attempt_num=2, angle_split=8,
)
print("CF₃ 替换后的原子数 / 键数 / conformer：", cf3_replaced.GetNumAtoms(), cf3_replaced.GetNumBonds(), cf3_replaced.GetNumConformers())
show_svg([parent, cf3_replaced], ["parent", "CF3 -> [*]COC(=O)C"])
show_3d([parent, cf3_replaced], ["parent", "CF3 -> [*]COC(=O)C"])

cf3_replaced = replace_substituent( parent, query="C(F)(F)F", replacement="[*]COC(=O)C", attempt_num=2, angle_split=8, ) print("CF₃ 替换后的原子数 / 键数 / conformer：", cf3_replaced.GetNumAtoms(), cf3_replaced.GetNumBonds(), cf3_replaced.GetNumConformers()) show_svg([parent, cf3_replaced], ["parent", "CF3 -> [*]COC(=O)C"]) show_3d([parent, cf3_replaced], ["parent", "CF3 -> [*]COC(=O)C"])

CF₃ 替换后的原子数 / 键数 / conformer： 107 113 1

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

3. bind_idx 消除多个匹配¶

c1ccccc1 在母体中有多个候选。下面选取连接到螯合碳的苯基，仅替换这一处。

In [4]:

phenyl_query = normalize_query("c1ccccc1")
phenyl_sites = find_sites(parent, phenyl_query)
carbon_bound_phenyl = next(
    site for site in phenyl_sites
    if parent.GetAtomWithIdx(site.scaffold_atom).GetSymbol() == "C"
)
bind_idx = carbon_bound_phenyl.scaffold_atom
bound_replaced = replace_substituent(
    parent, query="c1ccccc1", replacement="[*]CC#N",
    bind_idx=bind_idx, attempt_num=2, angle_split=8,
)
print("苯基候选数：", len(phenyl_sites), "；选中的 bind_idx：", bind_idx, "；conformer：", bound_replaced.GetNumConformers())
show_svg([parent, bound_replaced], ["parent", "selected phenyl -> [*]CC#N"])
show_3d([parent, bound_replaced], ["parent", "selected phenyl -> [*]CC#N"])

phenyl_query = normalize_query("c1ccccc1") phenyl_sites = find_sites(parent, phenyl_query) carbon_bound_phenyl = next( site for site in phenyl_sites if parent.GetAtomWithIdx(site.scaffold_atom).GetSymbol() == "C" ) bind_idx = carbon_bound_phenyl.scaffold_atom bound_replaced = replace_substituent( parent, query="c1ccccc1", replacement="[*]CC#N", bind_idx=bind_idx, attempt_num=2, angle_split=8, ) print("苯基候选数：", len(phenyl_sites), "；选中的 bind_idx：", bind_idx, "；conformer：", bound_replaced.GetNumConformers()) show_svg([parent, bound_replaced], ["parent", "selected phenyl -> [*]CC#N"]) show_3d([parent, bound_replaced], ["parent", "selected phenyl -> [*]CC#N"])

苯基候选数： 3 ；选中的 bind_idx： 5 ；conformer： 1

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

4. replace_all=True 批量替换¶

将所有符合悬挂片段契约的苯基替换为氰甲基；生成的原子会被内部标记，避免在循环中被重新匹配。

In [5]:

all_phenyl_replaced = replace_substituent(
    parent, query="c1ccccc1", replacement="[*]CC#N",
    replace_all=True, attempt_num=2, angle_split=8,
)
remaining_sites = find_sites(all_phenyl_replaced, phenyl_query)
print("批量替换后剩余符合条件的苯基位点：", len(remaining_sites), "；conformer：", all_phenyl_replaced.GetNumConformers())
show_svg([parent, all_phenyl_replaced], ["parent", "all pendant phenyl groups -> [*]CC#N"])
show_3d([parent, all_phenyl_replaced], ["parent", "all pendant phenyl groups -> [*]CC#N"])

all_phenyl_replaced = replace_substituent( parent, query="c1ccccc1", replacement="[*]CC#N", replace_all=True, attempt_num=2, angle_split=8, ) remaining_sites = find_sites(all_phenyl_replaced, phenyl_query) print("批量替换后剩余符合条件的苯基位点：", len(remaining_sites), "；conformer：", all_phenyl_replaced.GetNumConformers()) show_svg([parent, all_phenyl_replaced], ["parent", "all pendant phenyl groups -> [*]CC#N"]) show_3d([parent, all_phenyl_replaced], ["parent", "all pendant phenyl groups -> [*]CC#N"])

批量替换后剩余符合条件的苯基位点： 0 ；conformer： 1

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

5. 显式位点与金属配位键¶

令 query=None 并提供母体侧和取代基侧的原子索引，把真实的 Ni–C 配位键替换成 Ni–N 配位键。

In [6]:

nickel = next(atom for atom in parent.GetAtoms() if atom.GetSymbol() == "Ni")
nickel_carbon = next(neighbor for neighbor in nickel.GetNeighbors() if neighbor.GetSymbol() == "C")
dative_replaced = replace_substituent(
    parent, query=None, replacement="[*]N",
    start_idx=nickel.GetIdx(), end_idx=nickel_carbon.GetIdx(),
    attempt_num=2, angle_split=8,
)
print("显式替换的 Ni–C 原子索引：", nickel.GetIdx(), nickel_carbon.GetIdx(), "；conformer：", dative_replaced.GetNumConformers())
show_svg([parent, dative_replaced], ["parent", "explicit Ni-C -> Ni-N"])
show_3d([parent, dative_replaced], ["parent", "explicit Ni-C -> Ni-N"])

nickel = next(atom for atom in parent.GetAtoms() if atom.GetSymbol() == "Ni") nickel_carbon = next(neighbor for neighbor in nickel.GetNeighbors() if neighbor.GetSymbol() == "C") dative_replaced = replace_substituent( parent, query=None, replacement="[*]N", start_idx=nickel.GetIdx(), end_idx=nickel_carbon.GetIdx(), attempt_num=2, angle_split=8, ) print("显式替换的 Ni–C 原子索引：", nickel.GetIdx(), nickel_carbon.GetIdx(), "；conformer：", dative_replaced.GetNumConformers()) show_svg([parent, dative_replaced], ["parent", "explicit Ni-C -> Ni-N"]) show_3d([parent, dative_replaced], ["parent", "explicit Ni-C -> Ni-N"])

显式替换的 Ni–C 原子索引： 0 1 ；conformer： 1

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

6. 双键立体策略¶

Gaussian 母体没有适合作为孤立演示入口的非芳香双键，因此用最小的 FC=O 辅助分子展示 stereo_policy="E" 和 stereo_policy="Z"。

In [7]:

stereo_parent = Chem.MolFromSmiles("FC=O")
if stereo_parent is None:
    raise ValueError("无法构造双键立体演示分子。")
stereo_parent = Chem.AddHs(stereo_parent)
if AllChem.EmbedMolecule(stereo_parent, randomSeed=101) < 0:
    raise ValueError("无法为双键立体演示分子生成三维 conformer。")
stereo_e = replace_substituent(
    stereo_parent, query=None, replacement="[*]=[CH](Cl)",
    start_idx=1, end_idx=2, stereo_policy="E",
    attempt_num=2, randomSeed=101,
)
stereo_z = replace_substituent(
    stereo_parent, query=None, replacement="[*]=[CH](Cl)",
    start_idx=1, end_idx=2, stereo_policy="Z",
    attempt_num=2, randomSeed=101,
)
print("已生成 E 与 Z 两种连接双键构象。")
show_svg([stereo_parent, stereo_e, stereo_z], ["FC=O", "replacement | E", "replacement | Z"])
show_3d([stereo_parent, stereo_e, stereo_z], ["FC=O", "replacement | E", "replacement | Z"])

stereo_parent = Chem.MolFromSmiles("FC=O") if stereo_parent is None: raise ValueError("无法构造双键立体演示分子。") stereo_parent = Chem.AddHs(stereo_parent) if AllChem.EmbedMolecule(stereo_parent, randomSeed=101) < 0: raise ValueError("无法为双键立体演示分子生成三维 conformer。") stereo_e = replace_substituent( stereo_parent, query=None, replacement="[*]=[CH](Cl)", start_idx=1, end_idx=2, stereo_policy="E", attempt_num=2, randomSeed=101, ) stereo_z = replace_substituent( stereo_parent, query=None, replacement="[*]=[CH](Cl)", start_idx=1, end_idx=2, stereo_policy="Z", attempt_num=2, randomSeed=101, ) print("已生成 E 与 Z 两种连接双键构象。") show_svg([stereo_parent, stereo_e, stereo_z], ["FC=O", "replacement | E", "replacement | Z"]) show_3d([stereo_parent, stereo_e, stereo_z], ["FC=O", "replacement | E", "replacement | Z"])

已生成 E 与 Z 两种连接双键构象。

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

7. 实验性双锚点整体替换¶

多位点原型可以整体替换一个连通片段，也可以整体替换多个彼此独立的母体片段。这里用带三维构象的环氧样母体 C1CC2OC2C1，将 C1CO1 通过两个外部单键连接点替换为更大的并环片段。原型会冻结母体 scaffold 的坐标，只对新片段做刚体拟合和轴向旋转采样。

In [8]:

from molop.structure import replace_multisite_substituent
from molop.structure.substitution_sites import find_multisite_sites

multisite_parent = Chem.MolFromSmiles("C1CC2OC2C1")
if multisite_parent is None:
    raise ValueError("无法构造多锚点演示母体。")
multisite_parent = Chem.AddHs(multisite_parent)
if AllChem.EmbedMolecule(multisite_parent, randomSeed=7) < 0:
    raise ValueError("无法为多锚点演示母体生成三维 conformer。")

multisite_sites = find_multisite_sites(
    multisite_parent, normalize_query("C1CO1"), anchor_count=2
)
anchor_table = pd.DataFrame(
    [
        {
            "scaffold_idx": anchor.scaffold_atom,
            "old_fragment_idx": anchor.substituent_atom,
            "bond": str(anchor.bond_type),
        }
        for anchor in multisite_sites[0].anchors
    ]
)
display(anchor_table)
multisite_replaced = replace_multisite_substituent(
    multisite_parent,
    query="C1CO1",
    replacement="[*]C1C([*])C2CCC2C1",
    attempt_num=2,
    angle_split=36,
)
print("two-anchor replacement atoms / conformers:", multisite_replaced.GetNumAtoms(), multisite_replaced.GetNumConformers())
show_svg(
    [multisite_parent, multisite_replaced],
    ["epoxide-like parent", "larger fused ring"],
)
show_3d(
    [multisite_parent, multisite_replaced],
    ["epoxide-like parent", "larger fused ring"],
)

from molop.structure import replace_multisite_substituent from molop.structure.substitution_sites import find_multisite_sites multisite_parent = Chem.MolFromSmiles("C1CC2OC2C1") if multisite_parent is None: raise ValueError("无法构造多锚点演示母体。") multisite_parent = Chem.AddHs(multisite_parent) if AllChem.EmbedMolecule(multisite_parent, randomSeed=7) < 0: raise ValueError("无法为多锚点演示母体生成三维 conformer。") multisite_sites = find_multisite_sites( multisite_parent, normalize_query("C1CO1"), anchor_count=2 ) anchor_table = pd.DataFrame( [ { "scaffold_idx": anchor.scaffold_atom, "old_fragment_idx": anchor.substituent_atom, "bond": str(anchor.bond_type), } for anchor in multisite_sites[0].anchors ] ) display(anchor_table) multisite_replaced = replace_multisite_substituent( multisite_parent, query="C1CO1", replacement="[*]C1C([*])C2CCC2C1", attempt_num=2, angle_split=36, ) print("two-anchor replacement atoms / conformers:", multisite_replaced.GetNumAtoms(), multisite_replaced.GetNumConformers()) show_svg( [multisite_parent, multisite_replaced], ["epoxide-like parent", "larger fused ring"], ) show_3d( [multisite_parent, multisite_replaced], ["epoxide-like parent", "larger fused ring"], )

	scaffold_idx	old_fragment_idx	bond
0	1	2	SINGLE
1	5	4	SINGLE

two-anchor replacement atoms / conformers: 26 1

/tmp/ipykernel_2310/693395099.py:25: ExperimentalWarning: Experimental API 'replace_multisite_substituent' is not stable; its behavior and interface may change without notice.
  multisite_replaced = replace_multisite_substituent(

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

8. 三锚点与四锚点模式¶

三锚点和四锚点都来自同一个与母体骨架相连的环状片段：这里用一条连续的 scaffold 链分别连接目标三元环/四元环的全部原子，表格会逐行列出 3 条和 4 条边界键。此时原型不再绕单一连接轴采样，而是将全部锚点作为一个刚体几何整体拟合；替换侧和目标侧锚点都必须非共线。

In [9]:

def build_bridge_parent(anchor_count, random_seed):
    """Build a connected scaffold chain bridged to one multi-anchor ring."""
    rw_mol = Chem.RWMol()
    scaffold = [rw_mol.AddAtom(Chem.Atom("C")) for _ in range(anchor_count)]
    fragment = [
        rw_mol.AddAtom(Chem.Atom("N" if index == anchor_count - 1 else "C"))
        for index in range(anchor_count)
    ]
    for index in range(anchor_count - 1):
        rw_mol.AddBond(scaffold[index], scaffold[index + 1], Chem.BondType.SINGLE)
    for index in range(anchor_count):
        rw_mol.AddBond(fragment[index], fragment[(index + 1) % anchor_count], Chem.BondType.SINGLE)
        rw_mol.AddBond(scaffold[index], fragment[index], Chem.BondType.SINGLE)
    parent = rw_mol.GetMol()
    Chem.SanitizeMol(parent)
    parent = Chem.AddHs(parent)
    if AllChem.EmbedMolecule(parent, randomSeed=random_seed) < 0:
        raise ValueError("无法为多锚点演示母体生成三维 conformer。")
    return parent

three_anchor_parent = build_bridge_parent(anchor_count=3, random_seed=17)
four_anchor_parent = build_bridge_parent(anchor_count=4, random_seed=23)

three_anchor_sites = find_multisite_sites(
    three_anchor_parent, normalize_query("C1CN1"), anchor_count=3
)
four_anchor_sites = find_multisite_sites(
    four_anchor_parent, normalize_query("C1CCN1"), anchor_count=4
)
display(
    pd.DataFrame(
        [
            {
                "case": "three anchors",
                "anchor_rank": rank,
                "scaffold_idx": anchor.scaffold_atom,
                "fragment_idx": anchor.substituent_atom,
                "bond": str(anchor.bond_type),
            }
            for rank, anchor in enumerate(three_anchor_sites[0].anchors)
        ]
        + [
            {
                "case": "four anchors",
                "anchor_rank": rank,
                "scaffold_idx": anchor.scaffold_atom,
                "fragment_idx": anchor.substituent_atom,
                "bond": str(anchor.bond_type),
            }
            for rank, anchor in enumerate(four_anchor_sites[0].anchors)
        ]
    )
)
three_anchor_replaced = replace_multisite_substituent(
    three_anchor_parent,
    query="C1CN1",
    replacement="[*]C1C([*])C1([*])C2CC2",
    attempt_num=2,
    crowding_threshold=0.3,
)
four_anchor_replaced = replace_multisite_substituent(
    four_anchor_parent,
    query="C1CCN1",
    replacement="[*]C1C([*])C([*])C1([*])CCCC",
    attempt_num=2,
)
print("three-anchor / four-anchor atom counts:", three_anchor_replaced.GetNumAtoms(), four_anchor_replaced.GetNumAtoms())
show_svg(
    [three_anchor_parent, three_anchor_replaced, four_anchor_parent, four_anchor_replaced],
    ["3-anchor bridged parent", "3-anchor replacement", "4-anchor bridged parent", "4-anchor replacement"],
)
show_3d(
    [three_anchor_parent, three_anchor_replaced, four_anchor_parent, four_anchor_replaced],
    ["3-anchor bridged parent", "3-anchor replacement", "4-anchor bridged parent", "4-anchor replacement"],
)

def build_bridge_parent(anchor_count, random_seed): """Build a connected scaffold chain bridged to one multi-anchor ring.""" rw_mol = Chem.RWMol() scaffold = [rw_mol.AddAtom(Chem.Atom("C")) for _ in range(anchor_count)] fragment = [ rw_mol.AddAtom(Chem.Atom("N" if index == anchor_count - 1 else "C")) for index in range(anchor_count) ] for index in range(anchor_count - 1): rw_mol.AddBond(scaffold[index], scaffold[index + 1], Chem.BondType.SINGLE) for index in range(anchor_count): rw_mol.AddBond(fragment[index], fragment[(index + 1) % anchor_count], Chem.BondType.SINGLE) rw_mol.AddBond(scaffold[index], fragment[index], Chem.BondType.SINGLE) parent = rw_mol.GetMol() Chem.SanitizeMol(parent) parent = Chem.AddHs(parent) if AllChem.EmbedMolecule(parent, randomSeed=random_seed) < 0: raise ValueError("无法为多锚点演示母体生成三维 conformer。") return parent three_anchor_parent = build_bridge_parent(anchor_count=3, random_seed=17) four_anchor_parent = build_bridge_parent(anchor_count=4, random_seed=23) three_anchor_sites = find_multisite_sites( three_anchor_parent, normalize_query("C1CN1"), anchor_count=3 ) four_anchor_sites = find_multisite_sites( four_anchor_parent, normalize_query("C1CCN1"), anchor_count=4 ) display( pd.DataFrame( [ { "case": "three anchors", "anchor_rank": rank, "scaffold_idx": anchor.scaffold_atom, "fragment_idx": anchor.substituent_atom, "bond": str(anchor.bond_type), } for rank, anchor in enumerate(three_anchor_sites[0].anchors) ] + [ { "case": "four anchors", "anchor_rank": rank, "scaffold_idx": anchor.scaffold_atom, "fragment_idx": anchor.substituent_atom, "bond": str(anchor.bond_type), } for rank, anchor in enumerate(four_anchor_sites[0].anchors) ] ) ) three_anchor_replaced = replace_multisite_substituent( three_anchor_parent, query="C1CN1", replacement="[*]C1C([*])C1([*])C2CC2", attempt_num=2, crowding_threshold=0.3, ) four_anchor_replaced = replace_multisite_substituent( four_anchor_parent, query="C1CCN1", replacement="[*]C1C([*])C([*])C1([*])CCCC", attempt_num=2, ) print("three-anchor / four-anchor atom counts:", three_anchor_replaced.GetNumAtoms(), four_anchor_replaced.GetNumAtoms()) show_svg( [three_anchor_parent, three_anchor_replaced, four_anchor_parent, four_anchor_replaced], ["3-anchor bridged parent", "3-anchor replacement", "4-anchor bridged parent", "4-anchor replacement"], ) show_3d( [three_anchor_parent, three_anchor_replaced, four_anchor_parent, four_anchor_replaced], ["3-anchor bridged parent", "3-anchor replacement", "4-anchor bridged parent", "4-anchor replacement"], )

	case	anchor_rank	scaffold_idx	fragment_idx	bond
0	three anchors	0	0	3	SINGLE
1	three anchors	1	1	4	SINGLE
2	three anchors	2	2	5	SINGLE
3	four anchors	0	0	4	SINGLE
4	four anchors	1	1	5	SINGLE
5	four anchors	2	2	6	SINGLE
6	four anchors	3	3	7	SINGLE

three-anchor / four-anchor atom counts: 21 30

/tmp/ipykernel_2310/2746274293.py:54: ExperimentalWarning: Experimental API 'replace_multisite_substituent' is not stable; its behavior and interface may change without notice.
  three_anchor_replaced = replace_multisite_substituent(
/tmp/ipykernel_2310/2746274293.py:61: ExperimentalWarning: Experimental API 'replace_multisite_substituent' is not stable; its behavior and interface may change without notice.
  four_anchor_replaced = replace_multisite_substituent(

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

9. 断开的母体片段¶

多锚点 query 不要求自身连通。这里用 O.O 选中 OCCO 中两个彼此独立的氧端基，再用带两个标记的链整体替换它们；保留的母体骨架坐标不会移动。

In [10]:

disconnected_parent = Chem.MolFromSmiles("OCCO")
if disconnected_parent is None:
    raise ValueError("Could not construct the disconnected-query demo parent.")
disconnected_parent = Chem.AddHs(disconnected_parent)
if AllChem.EmbedMolecule(disconnected_parent, randomSeed=31) < 0:
    raise ValueError("Could not generate a 3D conformer for the disconnected-query demo.")

disconnected_sites = find_multisite_sites(
    disconnected_parent,
    normalize_query("O.O", require_connected=False),
    anchor_count=2,
)
if len(disconnected_sites) != 1:
    raise ValueError("The disconnected query did not resolve to one two-anchor site.")
disconnected_anchor_table = pd.DataFrame(
    [
        {
            "scaffold_idx": anchor.scaffold_atom,
            "old_fragment_idx": anchor.substituent_atom,
            "bond": str(anchor.bond_type),
        }
        for anchor in disconnected_sites[0].anchors
    ]
)
display(disconnected_anchor_table)
disconnected_replaced = replace_multisite_substituent(
    disconnected_parent,
    query="O.O",
    replacement="[*]CCCC[*]",
    attempt_num=2,
    crowding_threshold=0.1,
    anchor_tolerance=0.5,
    randomSeed=31,
)
print("Disconnected-query replacement atoms / conformers:", disconnected_replaced.GetNumAtoms(), disconnected_replaced.GetNumConformers())
show_svg(
    [disconnected_parent, disconnected_replaced],
    ["two independent parent fragments", "one bridged replacement"],
)
show_3d(
    [disconnected_parent, disconnected_replaced],
    ["two independent parent fragments", "one bridged replacement"],
)

disconnected_parent = Chem.MolFromSmiles("OCCO") if disconnected_parent is None: raise ValueError("Could not construct the disconnected-query demo parent.") disconnected_parent = Chem.AddHs(disconnected_parent) if AllChem.EmbedMolecule(disconnected_parent, randomSeed=31) < 0: raise ValueError("Could not generate a 3D conformer for the disconnected-query demo.") disconnected_sites = find_multisite_sites( disconnected_parent, normalize_query("O.O", require_connected=False), anchor_count=2, ) if len(disconnected_sites) != 1: raise ValueError("The disconnected query did not resolve to one two-anchor site.") disconnected_anchor_table = pd.DataFrame( [ { "scaffold_idx": anchor.scaffold_atom, "old_fragment_idx": anchor.substituent_atom, "bond": str(anchor.bond_type), } for anchor in disconnected_sites[0].anchors ] ) display(disconnected_anchor_table) disconnected_replaced = replace_multisite_substituent( disconnected_parent, query="O.O", replacement="[*]CCCC[*]", attempt_num=2, crowding_threshold=0.1, anchor_tolerance=0.5, randomSeed=31, ) print("Disconnected-query replacement atoms / conformers:", disconnected_replaced.GetNumAtoms(), disconnected_replaced.GetNumConformers()) show_svg( [disconnected_parent, disconnected_replaced], ["two independent parent fragments", "one bridged replacement"], ) show_3d( [disconnected_parent, disconnected_replaced], ["two independent parent fragments", "one bridged replacement"], )

	scaffold_idx	old_fragment_idx	bond
0	1	0	SINGLE
1	2	3	SINGLE

Disconnected-query replacement atoms / conformers: 18 1

/tmp/ipykernel_2310/232434793.py:26: ExperimentalWarning: Experimental API 'replace_multisite_substituent' is not stable; its behavior and interface may change without notice.
  disconnected_replaced = replace_multisite_substituent(

3Dmol.js failed to load for some reason. Please check your browser console for error messages.

10. 错误边界¶

任意子图替换并不是该 API 的目标。多匹配 query 必须使用 bind_idx 或 replace_all；环内键也会被明确拒绝。

In [11]:

errors = {}
try:
    replace_substituent(parent, query="C", replacement="[*]C", attempt_num=1)
except ValueError as error:
    errors["ambiguous query"] = str(error)
ring_bond = next(bond for bond in parent.GetBonds() if bond.IsInRing())
try:
    replace_substituent(
        parent, query=None, replacement="[*]C",
        start_idx=ring_bond.GetBeginAtomIdx(), end_idx=ring_bond.GetEndAtomIdx(),
        attempt_num=1,
    )
except ValueError as error:
    errors["ring bond"] = str(error)
display(pd.DataFrame(list(errors.items()), columns=["case", "error"]))

errors = {} try: replace_substituent(parent, query="C", replacement="[*]C", attempt_num=1) except ValueError as error: errors["ambiguous query"] = str(error) ring_bond = next(bond for bond in parent.GetBonds() if bond.IsInRing()) try: replace_substituent( parent, query=None, replacement="[*]C", start_idx=ring_bond.GetBeginAtomIdx(), end_idx=ring_bond.GetEndAtomIdx(), attempt_num=1, ) except ValueError as error: errors["ring bond"] = str(error) display(pd.DataFrame(list(errors.items()), columns=["case", "error"]))

	case	error
0	ambiguous query	query identifies multiple pendant substituents...
1	ring bond	substituent replacement does not support ring ...