AutoDock Vina +MGLtools (或AutoDockTools)是目前最主流的开源分子对接工具组合,由 Scripps 研究所开发,广泛应用于酶 - 底物结合预测、小分子虚拟筛选、蛋白 - 配体互作分析,完全适配各种合成生物学酶改造等研究场景,可用于评估突变对底物结合能力的影响、批量筛选优势突变体。本文主要介绍,如何从0开始在Linux系统中使用 AutoDock Vina 进行分子对接: | | | |
|---|
| AutoDock Vina | | 基于蒙特卡洛搜索 + 经验打分函数,预测小分子与蛋白的结合模式、结合亲和力 | |
| AutoDockTools (ADT) | | 将普通 PDB 格式的蛋白 / 配体转换为 Vina 专用的 PDBQT 格式(带电荷、原子类型、可旋转键信息);也用于对接盒子选取、结果可视化 | 有图形界面(GUI)和命令行脚本两种模式,服务器无桌面环境可用命令行脚本 |
适用场景包括但不限于:
conda create -n protein_simconda activate protein_sim
mamba install autodock-vina -ymamba install -c bioconda mgltools -ymamba install -c conda-forge openbabel=3.1.1 -y --override-channels
直接从 RCSB PDB 数据库下载原始结构(或者网页下载也可以):wget https://files.rcsb.org/download/9V6M.pdb
我们用 PET 降解的核心中间产物 MHET(单羟乙基对苯二甲酸)作为底物,SMILES 为:O=C(O)c1ccc(C(=O)OCCO)cc1
1. 蛋白受体预处理(9V6M.pdb → 9V6M_receptor.pdbqt)去除水分子、结晶杂原子,添加氢原子,计算 Gasteiger 电荷,转换为 Vina 专用的 PDBQT 格式(带电荷、原子类型信息)。prepare_receptor4.py \ -r 9V6M.pdb \ -o 9V6M_receptor.pdbqt \ -A hydrogens \ -U nphs_lps_waters
(1)9V6M 是单链 Apo 结构(无共晶配体),无需额外删配体,直接预处理即可(2)预处理后可以用以下命令检查原子数,确认蛋白完整grep ATOM 9V6M_receptor.pdbqt | wc -l
(3)若后续要做突变体对接,突变后的 PDB 文件用完全相同的参数预处理,保证变量唯一2. 小分子配体预处理(MHET → mhet.pdbqt)将 2D SMILES 转换为生理 pH 下的 3D 最优构象,添加氢、计算电荷,识别可旋转键,生成配体 PDBQT 文件。# 1. SMILES 转 3D 结构obabel -:"O=C(O)c1ccc(C(=O)OCCO)cc1" \ -O mhet.sdf \ --gen3d \ -p 7.4 # 生理pH 7.4下计算质子化状态# 2. SDF转MOL2格式obabel mhet.sdf -O mhet.mol2# 3. 用MOL2文件做配体预处理,生成 Vina 专用 PDBQTprepare_ligand4.py \ -l mhet.mol2 \ -o mhet.pdbqt \ -A hydrogens \ -U nphs_lps
对接 box 是 Vina 的搜索空间,必须完整覆盖活性位点,盒子偏移会直接导致结果无效。用 PyMOL 打开 9V6M.pdb,选中三个催化残基的 Cα 原子,读取坐标后取平均值:# 加载9V6M结构,替换为你本地的文件路径load 9V6M.pdb# 创建名为catalytic的选区,选中96/196/226号残基的Cα原子select catalytic, resi 96+196+226 and name CA# 在 PyMOL 命令行输入 Python 命令,直接输出平均坐标:pythonimport pymolcoords = pymol.cmd.get_coords("catalytic")center_x = (coords[0][0] + coords[1][0] + coords[2][0]) / 3center_y = (coords[0][1] + coords[1][1] + coords[2][1]) / 3center_z = (coords[0][2] + coords[1][2] + coords[2][2]) / 3print(f"对接Box中心:center_x = {center_x:.1f}, center_y = {center_y:.1f}, center_z = {center_z:.1f}")python end
resi 96+196+226:选择残基号为 96、196、226 的三个氨基酸(PETase 的催化三联体 Ser/Asp/His)name CA:只选中每个残基的 Cα 主链原子(三个原子刚好用来算中心)catalytic:自定义的选区名称,方便后续调用,执行后会看到蛋白上三个原子被高亮选中,确认位置在蛋白表面的活性口袋区域,说明残基号没选错。- 执行后会直接输出类似结果:对接盒子中心:center_x = -1.2, center_y = -0.4, center_z = -10.3,这三个数值就是你 Vina 配置文件里的盒子中心坐标。
MHET 是小分子底物,设置为20×20×20 埃即可(根据实际大小设置),完全覆盖活性口袋且不会引入过多无效搜索空间。
新建一个名为 vina_9v6m_config.txt 文件,所有参数对应先前建立的 9V6M+MHET 体系的参数:# 受体与配体文件receptor = 9V6M_receptor.pdbqtligand = mhet.pdbqt# 输出文件out = 9v6m_mhet_docking.pdbqtlog = 9v6m_mhet_log.txt# 对接盒子center_x = -1.2center_y = -0.4center_z = -10.3size_x = 20size_y = 20size_z = 20# 计算参数exhaustiveness = 32 # 保证结果稳定num_modes = 9 # 输出9个构象energy_range = 3 # 仅输出与最优构象能量差≤3kcal/mol的构象seed = 0 # 随机种子,保证结果可复现cpu = 16 # 调用CPU
vina --config vina_9v6m_config.txt
运行结束后终端会直接打印打分表,示例:
- 第 1 行是最优构象,结合能 - 5.5 kcal/mol,属于中等强度结合
结果解读:
1. 核心指标:结合能(Binding Affinity)- RMSD 表示当前构象与最优构象的结构差异(单位埃)
- 如果前 3 个构象的 RMSD 都 < 2 埃,说明对接结果稳定,置信度高。
3. 结果文件:
9v6m_mhet_docking.pdbqt:9 个对接构象的结构文件,可直接用 PyMOL 打开查看9v6m_mhet_log.txt:运行日志与打分表,可用于批量统计
# 加载受体和最优对接构象load 9V6M_receptor.pdbqt, receptorload 9v6m_mhet_docking.pdbqt, docking# 只显示第一个最优构象select mode1, docking and state 1# 查看催化残基与底物的距离select catalytic, receptor and resi 160+206+237dist catal_dist, catalytic, mode1
催化合理性验证:Ser 的羟基氧与底物酯键的距离在 3~5 埃范围内,才符合催化反应的空间要求。以上就是本次分享的所有内容,如果觉得有用,不妨点赞收藏转发给需要的人吧,这对我非常重要,非常感谢您的阅读!