smina-molecular-docking

smina 分子对接

概述

smina 是 AutoDock Vina 1.1.2 的一个分支，专注于灵活评分和能量最小化。可直接接受 SDF/MOL2/PDB 配体（无需手动转换为 PDBQT），可根据参考配体自动设置盒子，内置六种评分函数，并支持任意 --custom_scoring 项，还会输出逐原子评分贡献。仅提供 CLI——可通过 Python 的 subprocess 调用。

适用场景

• 对现有构象重新评分或进行局部能量最小化（--local_only、--minimize），无需完整搜索
• 不进行对接，仅计算单个构象的结合能（--score_only）
• 使用针对某类靶标调优的自定义或经验评分函数进行对接
• 输入 SDF/MOL2/多个配体，无需逐个将配体转换为 PDBQT
• 根据共结晶参考配体自动设置网格盒子
• 逐原子能量分解（--atom_term_data），用于药物化学类似物设计
• 适用于可在多个节点间高效并行的批量虚拟筛选（每个配体启动一个 CLI 进程）
• 如果需要 Vina Python 绑定、Vinardo 评分或 Vina 1.2 扩展后的力场，请改用 autodock-vina-docking；如果结合位点未知，请使用 diffdock

先决条件

• smina 二进制文件——conda-forge 或 SourceForge 构建版本
• Python：rdkit、openbabel-wheel（或系统 openbabel）、prody、pandas、py3Dmol
• ADFR Suite，用于 prepare_receptor（仅用于受体 PDBQT——配体由 smina 处理）
• 数据：蛋白质（PDB / PDB ID）、配体（SMILES / SDF / MOL2）

安装前请先检查——在 pixi/conda 环境中，smina 通常已在 PATH 中。如果 command -v smina 执行成功，则跳过安装；在 pixi 项目中，使用 pixi run smina ... 调用。

command -v smina || conda install -c conda-forge smina openbabel
pip install rdkit prody pandas py3Dmol
# ADFR Suite: https://ccsb.scripps.edu/adfr/downloads/

快速开始

使用参考配体自动设置盒子，进行端到端对接：

import subprocess

result = subprocess.run([
    "smina",
    "-r", "1hpv_receptor.pdbqt",
    "-l", "candidate.sdf",
    "--autobox_ligand", "1hpv_ref_ligand.pdb",
    "--autobox_add", "8",          # padding around reference (Å)
    "-o", "candidate_docked.sdf",
    "--exhaustiveness", "16",
    "--num_modes", "9",
    "--seed", "42",
], check=True, capture_output=True, text=True)

print(result.stdout.splitlines()[-15:])  # affinity table at stdout tail

工作流

步骤 1：准备受体（PDBQT）

移除水分子/杂原子，然后运行 ADFR Suite 的 prepare_receptor。

import subprocess, prody

pdb_id = "1HPV"
prody.fetchPDB(pdb_id, compressed=False)
protein = prody.parsePDB(f"{pdb_id}.pdb").select("protein")
prody.writePDB(f"{pdb_id}_protein.pdb", protein)

receptor_pdbqt = f"{pdb_id}_receptor.pdbqt"
subprocess.run([
    "prepare_receptor",
    "-r", f"{pdb_id}_protein.pdb",
    "-o", receptor_pdbqt,
    "-A", "hydrogens",
], check=True)
print(f"Receptor: {receptor_pdbqt} ({protein.numAtoms()} atoms)")

步骤 2：准备配体（SDF）

smina 可直接读取 SDF。使用 RDKit 生成 3D 坐标。

from rdkit import Chem
from rdkit.Chem import AllChem

mol = Chem.MolFromSmiles("CC(C)(C)NC(=O)[C@@H]1CN(CCc2ccccc2)C[C@H]1O")
mol = Chem.AddHs(mol)
AllChem.EmbedMolecule(mol, randomSeed=42)
AllChem.MMFFOptimizeMolecule(mol)

w = Chem.SDWriter("candidate.sdf"); w.write(mol); w.close()
print(f"Ligand SDF: candidate.sdf ({mol.GetNumAtoms()} atoms)")

步骤 3：提取参考配体以自动设置盒子

--autobox_ligand 根据参考结构推导网格。

import prody

ref = prody.parsePDB(f"{pdb_id}.pdb").select("hetero and not water and not ion")
if ref is None:
    raise RuntimeError("No reference ligand — supply explicit --center_x/--size_x")
prody.writePDB(f"{pdb_id}_ref_ligand.pdb", ref)
print(f"Ref ligand: {ref.numAtoms()} atoms, center {ref.getCoords().mean(axis=0).round(2)}")

步骤 4：使用自动设置的盒子运行对接

亲和力表会输出到 stdout——请将其捕获。

import subprocess

proc = subprocess.run([
    "smina",
    "-r", receptor_pdbqt,
    "-l", "candidate.sdf",
    "--autobox_ligand", f"{pdb_id}_ref_ligand.pdb",
    "--autobox_add", "8",
    "-o", "candidate_docked.sdf",
    "--exhaustiveness", "16",
    "--num_modes", "9",
    "--energy_range", "3",
    "--cpu", "4",
    "--seed", "42",
], check=True, capture_output=True, text=True)

for line in proc.stdout.splitlines()[-15:]:
    print(line)

步骤 5：解析构象和亲和力

亲和力会写入 SDF 的 <minimizedAffinity> 属性中。

from rdkit import Chem
import pandas as pd

rows = []
for i, mol in enumerate(Chem.SDMolSupplier("candidate_docked.sdf", removeHs=False)):
    if mol is None:
        continue
    aff = float(mol.GetProp("minimizedAffinity")) if mol.HasProp("minimizedAffinity") else None
    rmsd = float(mol.GetProp("minimizedRMSD")) if mol.HasProp("minimizedRMSD") else None
    rows.append({"pose": i + 1, "affinity_kcal_mol": aff, "rmsd_to_best": rmsd})

df = pd.DataFrame(rows).sort_values("affinity_kcal_mol")
print(df.to_string(index=False))
print(f"Best: {df.iloc[0]['affinity_kcal_mol']:.2f} kcal/mol")

步骤 6：对现有构象进行局部能量最小化

--local_only 可优化输入构象，无需全局搜索。

import subprocess

subprocess.run([
    "smina",
    "-r", receptor_pdbqt,
    "-l", "candidate_pose.sdf",
    "--autobox_ligand", "candidate_pose.sdf",   # box around the pose itself
    "--autobox_add", "4",
    "-o", "candidate_min.sdf",
    "--local_only",
    "--minimize_iters", "1000",
], check=True, capture_output=True, text=True)
print("Local minimization complete: candidate_min.sdf")

步骤 7：可视化对接复合物

import py3Dmol

with open(f"{pdb_id}_protein.pdb") as f: rec = f.read()
with open("candidate_docked.sdf") as f: lig = f.read()

view = py3Dmol.view(width=800, height=600)
view.addModel(rec, "pdb"); view.setStyle({"model": 0}, {"cartoon": {}})
view.addModel(lig.split("$$$$")[0] + "$$$$", "sdf")
view.setStyle({"model": 1}, {"stick": {}})
view.zoomTo({"model": 1}); view.show()

步骤 8：批量虚拟筛选

每个配体启动一个 smina 进程，可在多个核心或集群节点间高效并行。

import subprocess, pandas as pd
from rdkit import Chem
from rdkit.Chem import AllChem

library = pd.DataFrame({
    "name":   ["cpd_001", "cpd_002", "cpd_003"],
    "smiles": ["CC(=O)Oc1ccccc1C(=O)O",
               "CC(C)Cc1ccc(cc1)C(C)C(=O)O",
               "OC(=O)c1ccccc1O"],
})

results = []
for _, row in library.iterrows():
    lig_sdf, out_sdf = f"{row['name']}.sdf", f"{row['name']}_docked.sdf"
    mol = Chem.MolFromSmiles(row["smiles"]); mol = Chem.AddHs(mol)
    AllChem.EmbedMolecule(mol, randomSeed=42); AllChem.MMFFOptimizeMolecule(mol)
    w = Chem.SDWriter(lig_sdf); w.write(mol); w.close()

    proc = subprocess.run([
        "smina", "-r", receptor_pdbqt, "-l", lig_sdf,
        "--autobox_ligand", f"{pdb_id}_ref_ligand.pdb", "--autobox_add", "8",
        "-o", out_sdf, "--exhaustiveness", "8", "--num_modes", "1",
        "--cpu", "2", "--seed", "42",
    ], capture_output=True, text=True)

    if proc.returncode != 0:
        results.append({"name": row["name"], "affinity_kcal_mol": None})
        continue
    docked = next(iter(Chem.SDMolSupplier(out_sdf, removeHs=False)))
    aff = float(docked.GetProp("minimizedAffinity")) if docked and docked.HasProp("minimizedAffinity") else None
    results.append({"name": row["name"], "affinity_kcal_mol": aff})

ranked = pd.DataFrame(results).sort_values("affinity_kcal_mol")
ranked.to_csv("screening_results.csv", index=False)
print(ranked.to_string(index=False))

关键参数

参数	默认值	范围 / 选项	作用
--exhaustiveness	8	1-128	MC 搜索力度；生产环境使用 16-32，发表研究使用 64+
--num_modes	9	1-20	写入输出 SDF 的构象数
--energy_range	3	1-5 kcal/mol	输出中相对于最佳构象的最大 ΔE
--scoring	default	default、vina、vinardo、dkoes_scoring、dkoes_scoring_old、ad4_scoring	内置评分函数（默认值约等于 Vina）
--custom_scoring	—	路径	用户定义的加权项；会覆盖 --scoring
--autobox_ligand	—	PDB/SDF 路径	根据参考结构推导网格
--autobox_add	4	2-12 Å	自动设置的盒子范围周围的填充
--center_x/y/z、--size_x/y/z	—	Å	手动设置网格（自动设置盒子的替代方式）
--local_only	关闭	标志	跳过全局搜索；仅进行局部优化
--score_only	关闭	标志	仅计算能量，不进行能量最小化或搜索
--minimize	关闭	标志	进行能量最小化，不进行评分搜索
--minimize_iters	0（自动）	0-100000	局部优化器迭代次数
--atom_term_data	关闭	标志	在输出中提供逐原子评分贡献
--cpu	全部	1-N	每个任务使用的线程数
--seed	随机	整数	用于重现结果的 RNG 种子
--no_lig	关闭	标志	将仅含受体的评分作为基准

常用方案

配方：仅评分（单点能量）

无需重新对接即可比较不同评分函数下的构象能量。

import subprocess

proc = subprocess.run([
    "smina", "-r", "receptor.pdbqt", "-l", "pose.sdf",
    "--score_only", "--scoring", "vinardo",   # try vina, vinardo, dkoes_scoring
], capture_output=True, text=True, check=True)

for line in proc.stdout.splitlines():
    if line.startswith("Affinity:"):
        print(line)

配方：自定义评分函数

针对目标类别调优的经验评分；复现 dkoes 项或自定义拟合集。

import subprocess

with open("my_scoring.txt", "w") as f:
    f.write("""\
-0.035579    gauss(o=0,_w=0.5,_c=8)
-0.005156    gauss(o=3,_w=2,_c=8)
0.840245     repulsion(o=0,_c=8)
-0.035069    hydrophobic(g=0.5,_b=1.5,_c=8)
-0.587439    non_dir_h_bond(g=-0.7,_b=0,_c=8)
1.923        num_tors_div
""")

subprocess.run([
    "smina", "-r", "receptor.pdbqt", "-l", "candidate.sdf",
    "--custom_scoring", "my_scoring.txt",
    "--autobox_ligand", "ref.pdb", "--autobox_add", "8",
    "-o", "custom_docked.sdf", "--exhaustiveness", "16",
], check=True)

配方：逐原子评分分解

确定哪些配体原子对结合起主要作用。

import subprocess
from rdkit import Chem

subprocess.run([
    "smina", "-r", "receptor.pdbqt", "-l", "pose.sdf",
    "--score_only", "--atom_term_data", "-o", "pose_decomp.sdf",
], check=True)

mol = next(iter(Chem.SDMolSupplier("pose_decomp.sdf", removeHs=False)))
for prop in mol.GetPropNames():
    if "atom_term" in prop:
        print(f"{prop}: {mol.GetProp(prop)[:120]}")

配方：重新对接验证

重新对接共结晶配体；确认 RMSD < 2.0 Å。

import subprocess
from rdkit import Chem
from rdkit.Chem import AllChem

subprocess.run([
    "smina", "-r", "receptor.pdbqt", "-l", "ref_ligand.sdf",
    "--autobox_ligand", "ref_ligand.sdf", "--autobox_add", "8",
    "-o", "redocked.sdf", "--exhaustiveness", "32",
    "--num_modes", "1", "--seed", "42",
], check=True)

ref = next(iter(Chem.SDMolSupplier("ref_ligand.sdf", removeHs=False)))
docked = next(iter(Chem.SDMolSupplier("redocked.sdf", removeHs=False)))
rmsd = AllChem.GetBestRMS(ref, docked)
print(f"Re-docking RMSD: {rmsd:.2f} Å  ->  {'PASS' if rmsd < 2.0 else 'FAIL'}")

预期输出

• *_docked.sdf — 多模型 SDF，每个构象包含一个分子；亲和力位于 <minimizedAffinity> 中，至最佳构象的 RMSD 位于 <minimizedRMSD> 中
• screening_results.csv — 名称、SMILES、亲和力（kcal/mol）
• *_receptor.pdbqt — 已准备的受体
• 标准输出 — 亲和力表：构象、亲和力（kcal/mol）、RMSD 下限/上限

故障排除

问题	原因	解决方案
smina: command not found	二进制文件不在 PATH 中	conda install -c conda-forge smina；在 pixi 环境中使用 pixi run smina
Parse error on line ...（配体）	SDF/MOL2 格式错误	使用 RDKit SDWriter 在执行 AddHs + EmbedMolecule 后重新导出
Could not figure out box dimensions	缺少 --autobox_ligand 和 --center_x/--size_x	提供自动定框参考，或指定网格
亲和力为高度正值（>0）	配体位于网格外，或输入中存在空间冲突	增大 --autobox_add；确认参考构象位于结合位点内
--local_only 返回未更改的输入	在某些构建版本中，--minimize_iters 0 表示“不进行最小化”	显式设置 --minimize_iters 1000
多次运行得到相同构象	固定 --seed 且 --exhaustiveness 较低	提高到 32 以上；更改 --seed 以生成集成
受体 PDBQT 生成失败	prepare_receptor 不在 PATH 中	安装 ADFR Suite，并将 <adfr>/bin 添加到 PATH
对接后 SDF 为空	smina 静默中止	移除 capture_output，检查 proc.stderr；检查受体电荷
--score_only 与 --local_only 的评分不同	局部优化移动了构象	这是预期行为——score_only 用于固定构象，local_only 用于精修后的能量
拒绝 --custom_scoring 文件	项名称有误，或缺少权重列	每行格式：<weight><whitespace><term_name(args)> — 空白字符格式必须严格符合要求
逐原子分解结果为空	使用了 --minimize / 对接，而不是 --score_only	--atom_term_data 与 --score_only 搭配时最可靠

参考资料

• smina on SourceForge — 原始发行版
• smina GitHub mirror — 持续维护的构建版本
• Koes、Baumgartner、Camacho（2013） — smina + dkoes 评分论文，J Chem Inf Model
• Trott 与 Olson（2010） — 上游 Vina 论文，J Comput Chem
• ADFR Suite — 用于生成受体 PDBQT 的 prepare_receptor
• RDKit 文档 — 配体 3D 生成和 SDF I/O
• Open Babel — 另一种配体格式转换工具