深入解析蛋白质-小分子对接中的打分函数：从物理原理到AI智能

发表于2025-07-22|更新于2025-07-26|科技前沿

|浏览量:

你好，各位技术爱好者，我是你们的老朋友qmwneb946。今天，我们要潜入一个既古老又充满活力的交叉领域：计算生物学与药物发现。特别是，我们将聚焦于其中一个核心且最具挑战性的环节——蛋白质-小分子对接（Protein-Small Molecule Docking）中的“打分函数”（Scoring Function）。

想象一下，我们正在玩一个极其复杂的乐高积木游戏，其中一个巨大的积木代表着靶点蛋白，而无数微小的积木则代表着潜在的药物分子。我们的目标是找到那个能够完美契合蛋白“凹槽”的小分子，并与之紧密结合，从而发挥药效。这个寻找和评估“契合度”的过程，就是分子对接，而“契合度”的量化标准，便是打分函数。

在药物发现的漫长征程中，从浩瀚的化合物库中筛选出少数几个有希望的先导化合物，是极其耗时且昂贵的。计算方法，尤其是分子对接，极大地加速了这一过程，它让我们能够在计算机中预先筛选、优化甚至设计药物分子。而打分函数，正是分子对接的“灵魂”，它决定了我们能否准确地预测分子间的结合强度，并高效地识别出真正的“明星分子”。

这篇博客，我将带你深入剖析打分函数的各个方面。我们将从其物理化学基础谈起，逐步深入到各种类型（基于力场、经验性、知识库、以及当下最火热的机器学习/深度学习），探讨它们的原理、优缺点、应用场景，并展望未来的发展趋势。准备好了吗？让我们一起开启这段分子层面的探索之旅！

为什么我们需要打分函数？

在理解打分函数之前，我们首先要明白分子对接的整个流程。分子对接的核心任务是预测一个小分子（配体）与一个大分子（受体，通常是蛋白质）结合时的最佳姿态（构象和位置），并评估这种结合的强度。这个强度，通常用结合自由能（Binding Free Energy）来衡量。

分子对接算法首先会在受体结合口袋内尝试生成大量配体的不同构象和位置，这被称为“采样”或“搜索”过程。面对成千上万，甚至上亿种可能的结合姿态，我们如何判断哪种姿态是最稳定、最可能发生的？又如何比较不同配体与同一受体结合的倾向性？

这就是打分函数登场的时候了。打分函数是一个数学模型，它接受一个配体-受体复合物的构象作为输入，并输出一个数值。这个数值通常被解释为结合亲和力（affinity）或结合自由能的某种近似值。一个更低的（或有时是更高的，取决于约定）分数通常代表更强的结合亲和力。

打分函数在分子对接中扮演着双重角色：

姿态排名（Pose Ranking）：在搜索过程中，它用于评估并选择最佳的结合姿态。对于一个给定的配体，对接算法会生成多个可能的结合构象，打分函数会根据其分数对这些构象进行排序，选择分数最低（最优）的作为最终的预测结合姿态。
亲和力预测和虚拟筛选（Affinity Prediction & Virtual Screening）：对于多个不同的配体，打分函数用于预测它们与受体结合的相对亲和力，从而在大型化合物库中进行虚拟筛选，识别出最有可能成为药物的候选分子。

然而，构建一个准确、高效且普适的打分函数，是分子对接领域长期以来的一个巨大挑战。它不仅涉及到复杂的物理化学相互作用，还需应对生物大分子柔性、溶剂化效应、熵贡献等诸多复杂因素。

打分函数的分类：一场进化之旅

打分函数的研究和发展，是计算化学、生物信息学和机器学习等多学科交叉融合的成果。根据其背后的理论基础和构建方法，打分函数大致可以分为以下几类：

基于力场的打分函数

力场（Force Field）是描述原子间相互作用势能的数学模型。它们基于经典的物理学原理，如牛顿力学和电磁学，将分子系统中的原子视为带有电荷和质量的“球体”，并通过势能函数来描述它们之间的相互作用。基于力场的打分函数，试图直接计算配体-受体复合物的总能量，或其结合能。

基础原理：分子力场概述

一个典型的分子力场势能函数 $E_{\text{total}}$ 通常包括键合项和非键合项：

$E_{\text{total}} = E_{\text{bond}} + E_{\text{angle}} + E_{\text{torsion}} + E_{\text{non-bonded}}$

其中：

$E_{\text{bond}}$ ：键伸缩势能，描述共价键偏离平衡键长的能量惩罚。通常用简谐振子模型描述：
$E_{\text{bond}} = \sum_{\text{bonds}} k_b (r - r_0)^2$
其中 $k_b$ 是键力常数， $r$ 是当前键长， $r_0$ 是平衡键长。
$E_{\text{angle}}$ ：键角弯曲势能，描述键角偏离平衡键角的能量惩罚。同样常用简谐振子模型：
$E_{\text{angle}} = \sum_{\text{angles}} k_\theta (\theta - \theta_0)^2$
其中 $k_\theta$ 是键角力常数， $\theta$ 是当前键角， $\theta_0$ 是平衡键角。
$E_{\text{torsion}}$ ：二面角扭转势能，描述绕共价键旋转时的能量变化。通常用周期函数描述：
$E_{\text{torsion}} = \sum_{\text{dihedrals}} \frac{V_n}{2} [1 + \cos(n\phi - \delta)]$
其中 $V_n$ 是势垒高度， $n$ 是周期性， $phi$ 是二面角， $\delta$ 是相位角。
$E_{\text{non-bonded}}$ ：非键合相互作用势能，这是分子对接打分函数中最关键的部分。它主要包括范德华（van der Waals）相互作用和静电（Electrostatic）相互作用。
- 范德华相互作用（van der Waals）：包括色散力（诱导偶极-诱导偶极）、取向力（永久偶极-永久偶极）和诱导力（永久偶极-诱导偶极）。常用Lennard-Jones（LJ）12-6势能函数描述：
  $E_{\text{vdw}} = \sum_{i<j} \left[ A_{ij} \left(\frac{1}{R_{ij}}\right)^{12} - B_{ij} \left(\frac{1}{R_{ij}}\right)^{6} \right]$
  或更常见的形式：
  $E_{\text{vdw}} = \sum_{i<j} \epsilon_{ij} \left[ \left(\frac{R_{ij}^0}{R_{ij}}\right)^{12} - 2\left(\frac{R_{ij}^0}{R_{ij}}\right)^{6} \right]$
  其中 $R_{ij}$ 是原子 $i$ 和 $j$ 之间的距离， $A_{ij}$ 和 $B_{ij}$ 是参数， $\epsilon_{ij}$ 是势阱深度， $R_{ij}^0$ 是平衡距离。第一项代表短程排斥（电子云重叠），第二项代表长程吸引。
- 静电相互作用（Electrostatic）：描述带电原子之间的库仑相互作用。
  $E_{\text{electrostatic}} = \sum_{i<j} \frac{q_i q_j}{\epsilon R_{ij}}$
  其中 $q_i, q_j$ 是原子 $i, j$ 的部分电荷， $\epsilon$ 是介电常数， $R_{ij}$ 是原子间距离。

应用于打分函数

基于力场的打分函数通常计算配体-受体复合物的结合能 $\Delta E_{\text{binding}}$ ，它可以通过以下方式近似：

$\Delta E_{\text{binding}} = E_{\text{complex}} - (E_{\text{receptor}} + E_{\text{ligand}})$

其中 $E_{\text{complex}}$ 是配体-受体复合物的总能量， $E_{\text{receptor}}$ 是单独受体的能量， $E_{\text{ligand}}$ 是单独配体的能量。然而，这种简单的计算往往忽略了溶剂化效应和构象熵等重要贡献。

为了更准确地考虑溶剂化效应，常常会引入隐式溶剂模型（Implicit Solvation Models），例如广为人知的MM/PBSA（Molecular Mechanics / Poisson-Boltzmann Surface Area）和MM/GBSA（Molecular Mechanics / Generalized Born Surface Area）。这些方法将溶剂化自由能项加入到结合能计算中：

$\Delta G_{\text{binding}} = \Delta E_{\text{MM}} + \Delta G_{\text{solvation}} - T\Delta S_{\text{conformational}}$

其中：

$\Delta E_{\text{MM}}$ 是分子力学能量变化，包括键合和非键合项。
$\Delta G_{\text{solvation}}$ 是溶剂化自由能变化，通常分解为极性项（由Poisson-Boltzmann或Generalized Born方程计算）和非极性项（与溶剂可及表面积SA有关）。
$- T\Delta S_{\text{conformational}}$ 是构象熵贡献，这部分最难精确计算，在许多简化模型中常被忽略或粗略估计。

优点：

物理直观性强：基于原子间的实际物理相互作用，具有较好的可解释性。
普适性较好：参数相对稳定，对于不同类型的分子体系具有一定的通用性。
理论基础扎实：背后的物理模型经过了广泛验证。

缺点：

计算成本高昂：精确计算所有原子间的相互作用非常耗时，特别是对于大分子体系和大量原子。
参数化挑战：力场参数的准确性至关重要，但很难精确确定，尤其是一些特殊原子类型或相互作用。
忽略熵贡献：大多数简单的力场打分函数难以准确处理结合过程中的构象熵变化，这可能导致亲和力预测偏差。
溶剂化处理不足：显式溶剂模型计算量巨大，而隐式溶剂模型虽然高效但精度有限。
受限于构象采样：如果对接算法未能找到真正的最低能量构象，即使力场完美，打分也无法准确。

代表性软件/算法：

AMBER, CHARMM, OPLS：经典的分子力场本身。
AutoDock（一部分）
DOCK（早期版本）
Glide SP/XP：Schrödinger公司的对接软件Glide使用的打分函数，结合了经验性项和力场项，并加入了更多的物理修正。

一个简化的力场打分函数伪代码示例：

# 假设我们有一个包含原子信息（类型、坐标、电荷）的配体和受体对象

def calculate_energy(atom1, atom2, type_map):
    """
    计算两个原子之间的非键合相互作用能。
    这只是一个极度简化的示例，实际力场要复杂得多。
    """
    R_ij = calculate_distance(atom1.coords, atom2.coords)

    # Lennard-Jones (van der Waals)
    epsilon_ij, R0_ij = get_LJ_parameters(atom1.type, atom2.type, type_map)
    E_vdw = epsilon_ij * ((R0_ij / R_ij)**12 - 2 * (R0_ij / R_ij)**6)

    # Electrostatic (Coulomb)
    q1 = atom1.charge
    q2 = atom2.charge
    dielectric_constant = 1.0 # 真空或溶剂的介电常数，实际情况复杂
    E_electrostatic = (q1 * q2) / (dielectric_constant * R_ij)

    return E_vdw + E_electrostatic

def score_force_field(receptor_atoms, ligand_atoms, type_map):
    """
    基于力场计算配体-受体复合物的结合分数。
    这里的“分数”是结合能的近似值。
    """
    total_interaction_energy = 0.0

    # 遍历配体和受体中的所有原子对
    for lig_atom in ligand_atoms:
        for rec_atom in receptor_atoms:
            # 忽略同一分子内部的相互作用，这里只考虑配体-受体间
            # 并且需要处理1-2, 1-3, 1-4相互作用的特殊情况，这里简化
            total_interaction_energy += calculate_energy(lig_atom, rec_atom, type_map)

    # 实际中还需要考虑配体内部能量变化、溶剂化、熵等复杂因素
    # 返回负值代表更强的结合，所以结合能越低越好
    return total_interaction_energy

# type_map 示例 (实际是巨大的参数表)
# type_map = {
#     ('C_aromatic', 'O_carbonyl'): {'epsilon': 0.1, 'R0': 3.5, ...},
#     ...
# }

经验性打分函数

经验性打分函数（Empirical Scoring Functions）不像力场那样完全基于物理原理，而是通过对大量实验结合亲和力数据进行回归分析（线性或非线性），拟合出一个结合亲和力的经验公式。它们通常将结合能分解为一系列加和的化学项，每个项都对应一种特定的相互作用类型，并为其赋予一个拟合系数。

基础原理：线性回归模型

一个典型的经验性打分函数的形式通常如下：

$\Delta G_{\text{binding}} = \sum_i w_i f_i (\text{interaction})$

其中：

$\Delta G_{\text{binding}}$ 是预测的结合自由能。
$f_i$ 是第 $i$ 种相互作用的项，例如氢键、疏水效应、范德华力、静电相互作用、可旋转键的熵罚等。
$w_i$ 是对应项的权重或系数，通过回归分析（如最小二乘法）从实验数据中拟合得到。

常见的相互作用项包括：

氢键（Hydrogen Bonds, H-bonds）：通常根据氢键的几何学特性（距离、角度）来计算，可能区分供体和受体类型。
疏水相互作用（Hydrophobic Interactions）：通常与配体和受体之间接触的疏水表面积相关。
范德华相互作用（van der Waals）：通常是一个短程吸引和排斥的组合，可能比力场中的LJ势更简化。
静电相互作用（Electrostatic）：简化版库仑相互作用或仅考虑离子键。
配体柔性/熵（Ligand Flexibility/Entropy）：通常通过可旋转键的数量来惩罚，因为旋转键越多，结合时构象熵损失越大。
金属配位（Metal Coordination）：如果结合口袋中存在金属离子，会专门为此添加项。

示例表达式（伪）：
$\Delta G_{\text{binding}} = W_{\text{vdw}} \sum E_{\text{vdw}} + W_{\text{hbond}} \sum E_{\text{hbond}} + W_{\text{elec}} \sum E_{\text{elec}} + W_{\text{hydrophobic}} \times A_{\text{hydrophobic}} + W_{\text{rot}} \times N_{\text{rot}} + C$

其中 $N_{\text{rot}}$ 是配体中可旋转键的数量， $A_{\text{hydrophobic}}$ 是疏水接触面积， $C$ 是一个常数项。

构建流程

收集数据集：需要一个包含大量配体-受体复合物结构（PDB数据）以及对应的实验结合亲和力数据（如 $K_i, K_d, IC_{50}$ 等）的数据库。
特征提取：对于每个复合物，从结构中提取上述定义的各种相互作用项的数值。
回归拟合：使用回归算法（如多元线性回归、支持向量回归等）拟合各项的权重，使得预测值与实验值之间的误差最小。

优点：

计算速度快：各项计算相对简单，因此打分速度远快于基于力场的函数。
在训练数据上表现良好：如果训练数据具有代表性且质量高，在类似体系上通常能表现出良好的预测能力。
兼顾多种相互作用：通过不同项的权重，可以平衡各种化学相互作用的重要性。

缺点：

普适性有限：模型的准确性高度依赖于训练数据的多样性和质量。对于训练数据中未包含的相互作用类型或分子骨架，预测能力可能下降。
缺乏物理可解释性：权重是通过统计拟合得到的，不一定直接对应真实的物理常数，可能存在过拟合问题。
忽略协同效应：简单加和模型难以捕捉不同相互作用之间的复杂协同或拮抗效应。
不适用于从头设计：由于参数依赖于已知结合模式，难以用于预测全新的结合模式。

代表性软件/算法：

AutoDock Vina：虽然其名称包含“AutoDock”，但Vina的打分函数本质上是一个高度优化的经验性打分函数，结合了部分知识库的思想。
SCORE
FlexX
CDocker
ChemPLP

一个简化的经验性打分函数伪代码示例 (AutoDock Vina思想启发):

def calculate_vina_score(receptor_atoms, ligand_atoms, rotatable_bonds_count):
    """
    一个受AutoDock Vina启发的高度简化的经验性打分函数。
    Vina的真实打分函数更为复杂，包含多种项和加权。
    """
    score = 0.0

    # 1. 亲水/疏水相互作用 (这里简化为原子对的类型判断)
    # 假设我们有原子类型（C_hydrophobic, O_polar, N_polar等）
    hydrophobic_contact_count = 0
    polar_contact_energy = 0.0

    for lig_atom in ligand_atoms:
        for rec_atom in receptor_atoms:
            distance = calculate_distance(lig_atom.coords, rec_atom.coords)

            # 范德华/短程排斥：Vina使用高斯函数或线性衰减函数
            # 这里用一个简单的距离惩罚和吸引
            if distance < 2.0: # 短程排斥
                score += 10.0 * (2.0 - distance) # 距离越近惩罚越大
            elif 2.0 <= distance < 5.0: # 吸引范围
                score -= 0.5 * (5.0 - distance) # 距离越近吸引越大

            # 氢键：简化为满足距离和角度的极性原子对
            if is_polar_atom(lig_atom) and is_polar_atom(rec_atom) and 2.5 <= distance <= 3.5:
                # 实际的氢键项会考虑角度依赖性
                polar_contact_energy -= 5.0 # 假设一个常数贡献

            # 疏水相互作用：简化为疏水原子之间的接触
            if is_hydrophobic_atom(lig_atom) and is_hydrophobic_atom(rec_atom) and distance < 4.0:
                hydrophobic_contact_count += 1

    score += polar_contact_energy
    score -= 0.1 * hydrophobic_contact_count # 疏水接触越多，分数越低（越好）

    # 2. 可旋转键的熵罚
    score += 0.5 * rotatable_bonds_count # 旋转键越多，分数越高（越差）

    # 3. 经验常数项 (拟合得出)
    score += 10.0 # 调整分数基线

    # Vina的目标是找到最低分数，所以通常是负值，代表更强的亲和力
    return score

def is_polar_atom(atom):
    # 假设原子类型或电荷可以判断极性
    return atom.type in ['O', 'N', 'S'] or abs(atom.charge) > 0.2

def is_hydrophobic_atom(atom):
    # 假设原子类型可以判断疏水性
    return atom.type == 'C' and not is_polar_atom(atom)

知识库打分函数

知识库打分函数（Knowledge-Based Scoring Functions），也称为统计势能（Statistical Potentials），它们从大量已知蛋白质-配体复合物结构数据库（如PDB）中提取统计信息，推断出原子对之间在结合界面出现的概率分布。通过比较结合状态下原子对的出现频率与非结合状态（或随机状态）下的频率，来推导出一种“偏好势能”。

基础原理：反向玻尔兹曼法

其核心思想是，如果某种相互作用在已知的稳定复合物中频繁出现，那么这种相互作用就是有利的；反之，如果很少出现，则可能是不利的。这种偏好可以转化为“伪势能”或“统计势能”。

通常，一个知识库打分函数将结合能表示为一系列原子对势能之和：

$\Delta G_{\text{binding}} = \sum_{i \in \text{ligand}} \sum_{j \in \text{receptor}} E_{ij}^{\text{stat}}(d_{ij})$

其中 $E_{ij}^{\text{stat}}(d_{ij})$ 是原子类型 $i$ 和 $j$ 在距离 $d_{ij}$ 下的统计势能。

这个统计势能通常通过以下公式导出（基于玻尔兹曼分布的反推）：

$E_{ij}^{\text{stat}}(d) = -kT \ln \left( \frac{P_{ij}^{\text{bound}}(d)}{P_{ij}^{\text{reference}}(d)} \right)$

其中：

$k$ 是玻尔兹曼常数， $T$ 是温度。
$P_{ij}^{\text{bound}}(d)$ 是在真实蛋白质-配体复合物中，原子类型 $i$ 和 $j$ 在距离 $d$ 处出现的频率或概率密度。
$P_{ij}^{\text{reference}}(d)$ 是在随机或非结合状态下，原子类型 $i$ 和 $j$ 在距离 $d$ 处出现的参考频率或概率密度。参考态的构建至关重要，它通常通过随机抽样、非结合状态的原子对分布或均匀分布来模拟。

构建流程：

数据收集：从PDB等数据库中收集大量高质量的蛋白质-配体复合物结构。
原子类型定义：将原子分类为不同的类型（例如，根据元素、杂化态、连接的化学环境等），以捕捉更细致的相互作用差异。
距离统计：计算数据库中所有原子对在不同距离间隔（bin）内的出现频率，得到 $P_{ij}^{\text{bound}}(d)$ 。
参考态构建：构建一个参考数据集或模型来计算 $P_{ij}^{\text{reference}}(d)$ 。
势能导出：根据上述公式计算每个原子对在不同距离下的统计势能。

优点：

计算速度快：一旦势能表构建完成，打分过程只需查表和加和。
捕捉真实的生物相互作用：由于是基于真实生物结构的数据，它们能隐式地包含一些复杂的效应，如诱导契合、溶剂化、熵等，这些效应很难通过简单物理模型显式建模。
避免参数拟合：不需要进行复杂的回归拟合，直接从统计数据中推导。

缺点：

依赖于数据库质量和覆盖范围：如果数据库中某个特定的相互作用类型或原子对的例子很少，那么相应的统计势能将不准确。
假设静态结构：基于静态PDB结构，难以捕捉分子柔性和动态过程。
参考态的定义问题：如何准确定义“非结合”或“随机”状态是一个难题，不同的参考态会导致不同的势能。
忽略多体效应：通常只考虑两体相互作用，忽略了三个或更多原子同时参与的复杂相互作用。

代表性软件/算法：

PMF (Potential of Mean Force)
DrugScore
BLEEP (Bayesian LEarned Potentials)
SMoG (Statistical Mechanics of Gating)
ITScore

一个简化的知识库打分函数伪代码示例：

# 假设我们已经预先计算并存储了统计势能表
# statistical_potentials = {
#     ('atom_type_A', 'atom_type_B'): {
#         'bin_0-1A': energy_val_0,
#         'bin_1-2A': energy_val_1,
#         ...
#     },
#     ...
# }

def get_statistical_potential(atom_type1, atom_type2, distance, statistical_potentials):
    """
    从预计算的统计势能表中查找给定原子对和距离的势能。
    """
    # 确保原子类型顺序一致，例如按字母排序
    sorted_types = tuple(sorted((atom_type1, atom_type2)))

    if sorted_types not in statistical_potentials:
        return 0.0 # 未知原子对，假设无贡献或忽略

    # 确定距离所在的bin
    # 实际中会有更复杂的binning策略
    if distance < 1.0:
        distance_bin = 'bin_0-1A'
    elif 1.0 <= distance < 2.0:
        distance_bin = 'bin_1-2A'
    # ... 更多距离区间
    else:
        distance_bin = 'bin_far' # 超过一定距离假设无贡献

    return statistical_potentials[sorted_types].get(distance_bin, 0.0)

def score_knowledge_based(receptor_atoms, ligand_atoms, statistical_potentials):
    """
    基于知识库的打分函数，计算配体-受体复合物的总统计势能。
    """
    total_statistical_potential = 0.0

    for lig_atom in ligand_atoms:
        for rec_atom in receptor_atoms:
            distance = calculate_distance(lig_atom.coords, rec_atom.coords)
            # 获取原子类型（例如 'C_aromatic', 'O_carbonyl', 'N_donor'）
            lig_atom_type = get_atom_type(lig_atom)
            rec_atom_type = get_atom_type(rec_atom)

            potential = get_statistical_potential(lig_atom_type, rec_atom_type, distance, statistical_potentials)
            total_statistical_potential += potential

    return total_statistical_potential

机器学习/深度学习打分函数

近年来，随着大数据和人工智能技术的发展，机器学习（Machine Learning, ML）和深度学习（Deep Learning, DL）方法在分子对接打分函数领域展现出前所未有的潜力。它们不再依赖于预设的物理公式或统计分布，而是通过从海量的结构-亲和力数据中自动学习复杂的非线性关系。

基础原理：从特征工程到神经网络

ML/DL打分函数的核心是：

特征表示（Feature Representation）：如何将蛋白质-配体复合物的3D结构转化为机器学习模型可以理解的数值或向量。
模型构建（Model Construction）：选择合适的机器学习或深度学习模型来学习特征与结合亲和力之间的映射关系。

特征表示：

基于描述符/指纹（Descriptor/Fingerprint-based）：将配体和/或结合口袋的化学性质编码成向量。例如，基于拓扑的指纹、原子类型和相互作用类型的计数等。
基于网格（Grid-based）：将结合口袋划分为3D网格，并在每个网格点上编码原子类型、电荷、疏水性等属性。这类似于图像，可以直接输入卷积神经网络（CNN）。
基于图（Graph-based）：将分子表示为图，原子是节点，化学键是边。这能更好地捕捉分子的连接性和拓扑结构，特别适合图神经网络（GNN）。
基于相互作用指纹（Interaction Fingerprints）：专门编码配体和受体之间的相互作用类型（如氢键、盐桥、疏水接触等）。

机器学习模型：

经典ML算法：
- 支持向量机（SVM）
- 随机森林（Random Forest, RF）：通过构建大量决策树并集成它们的结果来提高预测精度。RF-Score就是一个著名的例子，它使用多种距离和相互作用特征作为输入。
- 梯度提升（Gradient Boosting, XGBoost, LightGBM）
  这些模型通常需要人工进行特征工程。
深度学习模型：
- 卷积神经网络（Convolutional Neural Networks, CNNs）：特别适用于处理网格状的3D结构数据，能够自动学习局部特征。例如，DeepDTA、DeepBind、GNINA（使用3D-CNN预测亲和力）等。
- 图神经网络（Graph Neural Networks, GNNs）：直接在分子图上操作，无需将3D结构转换为网格，能更好地捕捉原子间的非局部依赖关系。例如，D-score、GraphBAR等。
- 结合型模型（Hybrid Models）：结合多种架构和特征，如将CNN和GNN结合，或者将经典ML的特征与DL模型结合。

训练流程

大规模数据集：需要比传统经验性方法更大规模、更高质量的蛋白质-配体复合物结构和实验结合亲和力数据。ChEMBL、BindingDB、PDBbind等是常用数据集。
特征提取/表示：将每个复合物转换为模型所需的输入格式（描述符向量、3D网格、图结构等）。
模型训练：使用训练数据集对模型进行训练，通过最小化预测值与真实值之间的损失函数（如均方误差RMSE）来优化模型参数。
验证与测试：使用独立的验证集和测试集评估模型的泛化能力和预测精度。

优点：

高精度：能够学习复杂的非线性关系，在许多基准测试中超越传统方法。
自动特征学习：深度学习模型（特别是CNN和GNN）能够从原始数据中自动学习有用的特征，减少了人工特征工程的负担。
处理大数据：能够有效利用大规模的结构-亲和力数据。
捕捉复杂相互作用：理论上能够捕捉传统方法难以建模的复杂相互作用模式和协同效应。

缺点：

数据依赖性强：需要海量的、高质量的、多样化的训练数据。数据不足或偏差可能导致模型过拟合或泛化能力差。
缺乏可解释性：大多数ML/DL模型是“黑箱”，难以理解其内部决策过程，这在药物发现这种需要深入理解机制的领域是一个挑战。
计算资源需求大：模型训练通常需要强大的计算资源（GPU）。
泛化能力挑战：对于与训练数据分布差异很大的新靶点或新化学空间，模型的性能可能下降。
采样问题依然存在：ML/DL打分函数提高了打分的准确性，但它并不能解决分子对接算法在构象采样上的挑战。如果最优构象没有被采样到，再好的打分函数也无济于事。

代表性软件/算法：

RF-Score (Random Forest)
DEEPSCORE (CNN-based)
KDEEP (Multiple Kernels for Drug-Target Interaction)
Pafnucy (CNN-based)
GNINA (使用3D-CNN进行打分，并结合了新的搜索算法)
GraphBAR (GNN-based)
AffinityNet (GNN-based)

一个简化的ML/DL打分函数概念伪代码示例 (基于简单特征的Random Forest):

from sklearn.ensemble import RandomForestRegressor
import numpy as np

# 假设我们已经有了以下特征提取函数
def extract_features(receptor_atoms, ligand_atoms, rotatable_bonds_count):
    """
    从复合物结构中提取一系列数值特征。
    这些特征可以是基于物理化学的，也可以是统计的。
    """
    features = []

    # 1. 相互作用特征 (基于原子类型和距离)
    hbond_count = 0
    hydrophobic_contact_area = 0.0
    salt_bridge_count = 0
    # ... 更多相互作用类型

    for lig_atom in ligand_atoms:
        for rec_atom in receptor_atoms:
            distance = calculate_distance(lig_atom.coords, rec_atom.coords)

            # 示例：氢键计数
            if is_potential_hbond(lig_atom, rec_atom, distance):
                hbond_count += 1
            # 示例：疏水接触面积 (简化)
            if is_hydrophobic_contact(lig_atom, rec_atom, distance):
                hydrophobic_contact_area += 0.1 # 每次接触增加一个小的面积单位

    features.append(hbond_count)
    features.append(hydrophobic_contact_area)
    features.append(salt_bridge_count)

    # 2. 配体描述符
    features.append(len(ligand_atoms)) # 配体原子数
    features.append(rotatable_bonds_count) # 配体可旋转键数
    # ... 其他配体属性，如分子量、LogP等

    # 3. 结合口袋描述符
    # features.append(pocket_volume)
    # features.append(pocket_hydrophobicity)

    return np.array(features).reshape(1, -1) # 为模型准备好形状

# --- 训练阶段 (通常离线完成) ---
def train_ml_scoring_function(training_data):
    """
    训练一个机器学习打分函数。
    training_data: 列表，每个元素是 (receptor_atoms, ligand_atoms, rotatable_bonds_count, experimental_affinity)
    """
    X_train = [] # 特征矩阵
    y_train = [] # 实验亲和力

    for rec_atoms, lig_atoms, rot_bonds, affinity in training_data:
        features = extract_features(rec_atoms, lig_atoms, rot_bonds)
        X_train.append(features[0]) # features是二维数组，取第一行
        y_train.append(affinity)

    # 使用Random Forest Regressor模型
    model = RandomForestRegressor(n_estimators=100, random_state=42, n_jobs=-1)
    model.fit(X_train, y_train)

    return model

# --- 预测阶段 ---
def predict_with_ml_model(model, receptor_atoms, ligand_atoms, rotatable_bonds_count):
    """
    使用训练好的ML模型预测结合亲和力。
    """
    features = extract_features(receptor_atoms, ligand_atoms, rotatable_bonds_count)
    predicted_affinity = model.predict(features)[0] # predict返回数组

    return predicted_affinity

# 辅助函数 (简化)
def calculate_distance(coords1, coords2):
    return np.linalg.norm(np.array(coords1) - np.array(coords2))

def is_potential_hbond(atom1, atom2, distance):
    # 极度简化：检查是否是H-donor/acceptor，且距离在合理范围内
    return (atom1.is_donor and atom2.is_acceptor or atom1.is_acceptor and atom2.is_donor) and 2.5 <= distance <= 3.5

def is_hydrophobic_contact(atom1, atom2, distance):
    # 极度简化：检查是否是疏水原子，且距离在合理范围内
    return atom1.is_hydrophobic and atom2.is_hydrophobic and distance <= 4.0

混合打分函数

许多现代的分子对接软件和打分函数并不是纯粹的某种类型，而是结合了多种方法的优点，形成了混合打分函数（Hybrid Scoring Functions）。例如，它们可能在经验性框架中融入力场计算的范德华和静电项，同时通过机器学习优化这些项的权重或引入额外的统计项。

这种混合方法旨在取长补短：利用力场的物理基础、经验方法的快速性和知识库的统计洞察力，并通过机器学习来更智能地学习和组合这些信息。

打分函数的挑战与局限

尽管打分函数取得了显著进步，但它们仍然面临诸多挑战：

准确性与普适性的权衡：
- 姿态预测（Pose Prediction）：打分函数需要准确识别出配体的正确结合构象。这通常比亲和力预测更容易，因为通常有多个“近乎正确”的构象。
- 亲和力预测（Affinity Prediction）：更具挑战性，需要精确量化结合自由能。目前的打分函数在预测精确结合能（例如，在 $\pm 1 \text{ kcal/mol}$ 的范围内）方面仍有很大差距。
- 排名能力（Ranking Capability）：对于不同分子，能否准确预测其相对亲和力并进行排序，这对于虚拟筛选至关重要。
蛋白质柔性与诱导契合：大多数对接算法和打分函数假定蛋白质是刚性的，或只允许局部柔性。然而，在配体结合时，蛋白质通常会发生构象变化（诱导契合），这会显著影响结合亲和力。精确建模这种动态过程是计算量巨大的挑战。
溶剂化效应：水分子在结合界面的作用至关重要。它们可以形成氢键、介导相互作用，甚至被排斥出结合口袋。准确计算脱溶剂化能（desolvation energy）和溶剂分子的熵贡献是一个复杂问题。
熵贡献的评估：结合过程不仅涉及能量变化，还涉及熵变化。配体和蛋白质在结合时自由度减少，导致构象熵损失。水分子被排斥出结合口袋时，其自由度增加，贡献正熵。准确量化这些熵项是极其困难的，特别是构象熵。
非加和性效应：打分函数通常将总结合能分解为独立相互作用项的简单加和。然而，实际的分子相互作用可能存在协同或拮抗效应，导致非加和性。
数据质量与偏差：对于经验性、知识库和机器学习打分函数，训练数据的质量、数量和多样性直接决定了模型的性能。数据中的误差、偏差或缺乏代表性都可能导致模型表现不佳。
“黑箱”问题：特别是对于复杂的ML/DL模型，其内部决策机制不透明，使得研究人员难以理解为何某些分子得分高或低，从而限制了药物设计中的迭代优化。

打分函数的发展趋势与未来展望

尽管面临挑战，打分函数的研究仍在飞速发展。未来的趋势主要集中在以下几个方面：

更强大的机器学习和深度学习模型：
- 更精细的特征表示：探索更有效的方式来编码3D结构和相互作用，例如，结合物理信息和学习特征。
- 图神经网络的深化：GNN在分子表示方面具有天然优势，未来将有更多结合物理原理和生物背景的GNN架构出现。
- 生成模型（Generative Models）：结合生成对抗网络（GANs）或变分自编码器（VAEs）等，不仅预测亲和力，还能生成具有高亲和力的分子。
- 多任务学习（Multi-task Learning）：同时预测结合亲和力、ADMET性质（吸收、分布、代谢、排泄、毒性）等，利用不同任务之间的相关性提高预测精度。
更准确的柔性处理：
- 集成分子动力学（MD）模拟：利用MD模拟来探索蛋白质和配体的构象空间，并从中提取更准确的相互作用模式和熵贡献。例如，通过MM/GBSA或MM/PBSA与MD结合。
- 更高效的柔性对接算法：开发能处理更大范围蛋白柔性的对接算法，例如，基于集成（ensemble docking）或诱导契合（induced fit docking）的方法。
- 可变形网格/图表示：在ML/DL模型中引入能够感知或适应分子形变的表示。
更好的溶剂化模型：
- 开发更精确且计算高效的隐式溶剂模型。
- 在ML/DL模型中更显式地考虑水分子在结合界面的作用。
可解释人工智能（XAI）：
- 为“黑箱”ML/DL模型开发可解释性工具，例如，通过显著性图（saliency maps）识别对预测贡献最大的原子或相互作用。这将有助于药物化学家理解分子结合机制，并指导分子优化。
数据驱动与物理模型的融合：
- 结合经典力场和量子化学计算的物理原理，与机器学习的数据学习能力。例如，使用ML来修正力场参数，或学习力场中难以建模的项。
- 通过神经网络势（Neural Network Potentials）直接从量子力学数据中学习原子间相互作用。
大规模、高质量的数据集建设：
- 持续扩充和完善蛋白质-配体结构和亲和力数据库，特别是那些具有高精度结构和多样化化学空间的复合体。
- 开发更智能的数据筛选和预处理方法，减少数据噪声。
集成与自动化：
- 将打分函数与其他药物发现工具（如分子动力学、量子化学、合成路线预测）更紧密地集成，形成端到端的自动化工作流。

结论

蛋白质-小分子对接中的打分函数，是计算药物发现领域的基石，也是连接理论物理化学与实际生物问题的桥梁。从早期基于物理力场的尝试，到经验性模型的快速发展，再到如今人工智能浪潮下机器学习和深度学习的崛起，打分函数一直在不断进化，以期更准确、更高效地预测分子间的“契合度”。

我们看到，没有哪一种打分函数是完美的。力场打分函数提供了物理基础但计算昂贵且难处理熵；经验性函数快速有效但泛化能力受限；知识库函数捕捉统计规律但依赖数据质量；而ML/DL虽然精度喜人，却面临数据依赖、可解释性差和计算资源高的挑战。

然而，正是这些挑战推动着科研人员不断探索和创新。未来的打分函数将不再局限于单一范式，而是走向多模型融合、数据驱动与物理原理深度结合的道路。随着计算资源的进步、更大规模高质量数据的积累，以及更智能算法的涌现，我们有理由相信，打分函数将变得越来越强大，最终成为药物发现和设计中不可或缺的“智慧之眼”，加速更多创新药物的诞生，为人类健康福祉贡献更大的力量。

希望这篇长文能让你对打分函数有了更深入的理解。如果你有任何问题或想法，欢迎在评论区与我交流！我是qmwneb946，我们下次再见！

文章作者: qmwneb946

文章链接: https://qmwneb946.dpdns.org/2025/07/22/2025-07-22-232848/

科技前沿 2025 蛋白质-小分子对接的打分函数