R-group Decomposition：从取代基拆解理解 SAR

R 基团分解把一组同骨架的分子拆成「共同核心 + 各自的取代基」，得到一张 R 基团表。这是药物化学分析 SAR 最直接的工具——它把「一堆分子和活性数据」变成「哪个位置换什么基团更好」。

为什么这个视角重要

# 【原始数据的形式】：
#   SMILES                          pIC50
#   CC(=O)Nc1ccc(F)cc1              6.2
#   CC(=O)Nc1ccc(Cl)cc1             6.8
#   CC(=O)Nc1ccc(Br)cc1             7.1
#   CC(=O)Nc1ccc(OC)cc1             5.9
#   ...
#   → 【看不出规律】
#
# 【R 基团分解后】：
#   核心: CC(=O)Nc1ccc([*:1])cc1
#   R1     pIC50
#   F      6.2
#   Cl     6.8
#   Br     7.1
#   OMe    5.9
#   → 【规律立刻显现】：
#     卤素随体积增大活性上升
#     甲氧基不如卤素
#     → 提示该位置需要疏水的、体积适中的基团
#
# 【这就是 SAR 分析的核心工作方式】（见 073）

用 RDKit 做 R 基团分解

from rdkit import Chem
from rdkit.Chem import rdRGroupDecomposition as rgd
import pandas as pd

# ---- 方法一：指定核心 ----
core = Chem.MolFromSmarts("CC(=O)Nc1ccccc1")
mols = [Chem.MolFromSmiles(s) for s in smiles_list]
mols = [m for m in mols if m is not None]

params = rgd.RGroupDecompositionParameters()
params.removeHydrogensPostMatch = True
params.onlyMatchAtRGroups = False

decomp = rgd.RGroupDecomposition(core, params)
matched_idx = []
for i, m in enumerate(mols):
    if decomp.Add(m) >= 0:
        matched_idx.append(i)
decomp.Process()

result = decomp.GetRGroupsAsColumns(asSmiles=True)
df = pd.DataFrame(result)
print(df.head())
# 列：Core, R1, R2, ...

print(f"匹配了 {len(matched_idx)}/{len(mols)} 个分子")
# 【未匹配的分子需要单独看】——可能属于其它系列

# ---- 方法二：自动找核心（用 MCS）----
from rdkit.Chem import rdFMCS

mcs = rdFMCS.FindMCS(
    mols,
    ringMatchesRingOnly=True,        # 【环只匹配环】
    completeRingsOnly=True,          # 【只匹配完整的环】
    timeout=60,
)
core = Chem.MolFromSmarts(mcs.smartsString)
print("最大公共子结构:", mcs.smartsString)

# 【MCS 参数的影响很大】：
#   completeRingsOnly=False 可能得到「半个环」的核心
#   → 化学上不合理
#   → 【做 SAR 分析时应该设为 True】

# ---- 方法三：用骨架作为核心 ----
from rdkit.Chem.Scaffolds import MurckoScaffold
core = MurckoScaffold.GetScaffoldForMol(mols[0])
# 【前提是这批分子确实同骨架】（见 041）

SAR 分析的实用流程

import pandas as pd
import numpy as np

# 把 R 基团表与活性数据合并
df_rgroup = pd.DataFrame(result)
df_rgroup["pIC50"] = [activities[i] for i in matched_idx]

# ---- 分析一：单个 R 位置的影响 ----
for col in [c for c in df_rgroup.columns if c.startswith("R")]:
    stats = df_rgroup.groupby(col)["pIC50"].agg(["count", "mean", "max", "std"])
    stats = stats[stats["count"] >= 2].sort_values("mean", ascending=False)
    print(f"\n=== {col} ===")
    print(stats.head(10))

# 【判读】：
#   某个位置的取代基之间活性差异大 → 【该位置敏感，是优化的关键】
#   差异小 → 该位置不敏感，可以用来调节 ADMET 性质
#   → 【这个区分对优化策略很重要】：
#     在不敏感的位置改 logP，不会损失活性

# ---- 分析二：找出「组合效应」----
#   两个 R 位置的组合是否有协同或拮抗？
pivot = df_rgroup.pivot_table(
    index="R1", columns="R2", values="pIC50", aggfunc="mean")
print(pivot)
# 【如果某个 R1-R2 组合明显好于各自的边际效应之和
#   → 存在协同，可能是特定的相互作用】

# ---- 分析三：找活性悬崖 ----
#   结构相似但活性差异大的取代基对
def find_cliffs(df, r_col, act_col="pIC50", threshold=1.0):
    from rdkit import Chem, DataStructs
    from rdkit.Chem import rdFingerprintGenerator
    gen = rdFingerprintGenerator.GetMorganGenerator(radius=2, fpSize=2048)

    groups = df.groupby(r_col)[act_col].mean().to_dict()
    items = list(groups.items())
    cliffs = []
    for i in range(len(items)):
        for j in range(i + 1, len(items)):
            r1, a1 = items[i]
            r2, a2 = items[j]
            if abs(a1 - a2) < threshold:
                continue
            m1 = Chem.MolFromSmiles(r1.replace("[*:1]", "C"))
            m2 = Chem.MolFromSmiles(r2.replace("[*:1]", "C"))
            if m1 is None or m2 is None:
                continue
            sim = DataStructs.TanimotoSimilarity(
                gen.GetFingerprint(m1), gen.GetFingerprint(m2))
            if sim > 0.6:
                cliffs.append((r1, r2, a1, a2, sim))
    return sorted(cliffs, key=lambda x: -abs(x[2] - x[3]))

# 【活性悬崖是信息量最大的数据点】：
#   它们指向具体的相互作用
#   → 应该重点分析并用结构解释（见 096）

可视化 R 基团表

from rdkit.Chem import Draw
from rdkit.Chem.Draw import rdMolDraw2D

def draw_rgroup_series(df, r_col, act_col="pIC50", top_n=12):
    """按活性排序，画出取代基"""
    stats = df.groupby(r_col)[act_col].mean().sort_values(ascending=False)
    mols, legends = [], []
    for r, act in list(stats.items())[:top_n]:
        # R 基团的 SMILES 含 [*:1] 标记
        m = Chem.MolFromSmiles(r)
        if m is None:
            continue
        mols.append(m)
        legends.append(f"{act:.2f}")
    return Draw.MolsToGridImage(mols, molsPerRow=4, legends=legends,
                                subImgSize=(220, 180))

# 【可视化的价值】：
#   化学家看图比看表格快得多
#   → 能立刻看出「大的疏水基团好」这类规律
#
# 【更好的工具】：
#   mols2grid（见 217）：交互式的分子表格
#   → 可以排序、筛选、加注释
#   → 【在 Notebook 中做 SAR 分析很方便】

实践中的问题

• 核心的选择影响很大：核心太大 → 匹配的分子少；核心太小 → R 基团太复杂难以解读。通常需要试几个；
• 对称性问题：如果核心有对称性（如对位取代的苯环），R1 与 R2 的分配可能不唯一——RDKit 会尝试处理，但需要检查；
• 未匹配的分子：它们可能属于其它系列，应该单独分析而非丢弃；
• 数据量要求：每个 R 基团至少要有 2～3 个数据点才有统计意义；
• 混淆因素：如果 R1 与 R2 同时变化，单独看 R1 的效应可能被 R2 干扰——理想的 SAR 数据是「一次只变一个位置」；
• 与 MMPA 的关系：MMPA（见 043《MMPA 配对分子分析》）是更自动化的方法，不需要预先指定核心。

延伸资源

• 骨架：040《Scaffold 骨架》、041《Bemis–Murcko Scaffold》；MMPA：043《MMPA 配对分子分析》；SAR：073《SAR 构效关系》；
• 可视化：217《Mols2grid》；姿势验证：096《Binding Pose》。