简介这份资料面向生物信息学与计算化学方向的研究者、研究生及药物设计初学者聚焦于在Ubuntu 18.04环境下搭建Autodock Vina分子对接流程并结合Slurm调度器实现集群上的批量对接任务。内容涵盖依赖库安装、Open Babel与MGLTools配置、受体与配体预处理脚本的使用以及通过conf.txt控制CPU数量与搜索空间等关键参数最终输出pdbqt并转换为sdf便于后续分析同时详细说明Slurm的munge认证、slurmctld与slurmd服务启动及slurm.conf编写方法。资源包为1个PDF文件大小约1013KB以图文步骤形式呈现完整操作链路。目前已有5185人学习适合希望从单机对接进阶到集群批量计算、提升虚拟筛选效率的读者参考。1. 从一次 200 个配体的对接任务说起单个分子对接跑一次 Autodock Vina 只要几十秒很多人第一次接触时觉得这东西挺快。但真正落到虚拟筛选场景比如你手上有 200 个从天然产物库筛出来的候选配体要对着同一个靶点蛋白逐个跑手动改配置、改文件名、点运行、等结果一天就搭进去了。批量分子对接要解决的核心问题不是怎么对接而是怎么把重复劳动交给脚本同时保证每个配体的参数一致、结果可追溯。Autodock Vina 本身是命令行工具这个特性决定了它天然适合被批处理脚本调用。它接受--receptor、--ligand、--center_x/y/z、--size_x/y/z这些参数输出一个包含若干结合构象和亲和力的文本结果。批量化的关键就在于把配体列表、受体文件、对接盒子参数抽成配置用一个循环或并行框架去驱动 Vina再把散落的输出统一解析成表格。适合读这篇的人是已经会跑单次对接、现在要把规模从 1 提到 100 以上的生信或药物设计从业者。2. Autodock Vina 批量对接前的文件准备与参数确定2.1 受体与配体的格式转换从 PDB 到 PDBQTVina 不直接吃 PDB 文件受体和配体都得转成 PDBQT。受体转换用 MGLTools 里的prepare_receptor4.py配体用prepare_ligand4.py。这一步是批量流程里最容易翻车的地方因为配体来源五花八门有的带氢有的不带有的电荷没加。# 受体准备去水、加氢、加电荷输出 PDBQT prepare_receptor4.py -r receptor.pdb -o receptor.pdbqt \ -A hydrogens -U nphs_lps_waters # 单个配体准备 prepare_ligand4.py -l ligand.sdf -o ligand.pdbqt -A hydrogens-A hydrogens表示自动加氢-U nphs_lps_waters表示删除非极性氢、孤立水分子和孤对电子。批量场景下配体往往是 SDF 或 MOL2 格式需要先拆分再逐个转换。常见做法是用 Open Babel 先把多分子 SDF 拆成单文件# 把含多个分子的 SDF 拆成独立文件 obabel library.sdf -O lig_.mol2 -m-m是 split 模式会生成lig_1.mol2、lig_2.mol2这样的序列文件。拆完之后再写个循环批量转 PDBQT比一个个手点靠谱得多。2.2 对接盒子中心与尺寸怎么定才不白跑对接盒子grid box决定了搜索空间。盒子太小配体可能被截断在边缘构象不合理盒子太大搜索空间爆炸Vina 跑得慢还容易给出假阳性。定盒子的标准做法是围绕已知活性位点或共晶配体中心。如果你有共晶结构最省事的办法是算共晶配体的几何中心# 用 MDAnalysis 算共晶配体中心作为盒子中心 import MDAnalysis as mda u mda.Universe(complex.pdb) lig u.select_atoms(resname LIG) # 换成你的配体残基名 center lig.center_of_mass() print(center_x , round(center[0], 3)) print(center_y , round(center[1], 3)) print(center_z , round(center[2], 3))center_of_mass()返回质心坐标直接填进 Vina 的--center_x/y/z。盒子尺寸一般让每个方向比配体最大跨度多留 8–12 Å常见取值是 20×20×20 到 30×30×30。没有共晶结构时用 fpocket、SiteMap 之类的工具预测口袋或者干脆把整个蛋白包进去做盲对接但盲对接的假阳性率会明显上升。2.3 参数表批量对接必须统一的几个值批量对接最忌讳每个配体用不同参数结果没法横向比较。下面这张表是我一般会固定下来的参数写进配置文件里让所有配体共用。参数含义常用取值说明--exhaustiveness搜索彻底程度8–32批量筛选用 8精确重跑用 32--num_modes输出构象数9默认 9够看排名--energy_range构象能量差上限3–4单位 kcal/mol--cpu使用核数按机器定单任务别占满留给并行--seed随机种子固定值保证可复现--exhaustiveness是批量场景里最需要权衡的参数。它和运行时间近似线性8 和 32 的耗时差 4 倍左右但打分排名在前几位的构象往往差别不大。做初筛时用 8 快速过一遍挑出 top 10% 再用 32 精算是性价比比较高的策略。3. 用 Python 驱动 Vina 实现批量分子对接3.1 单次 Vina 命令的拆解先把单次命令写清楚批量才有模板可套。一条完整的 Vina 命令长这样vina --receptor receptor.pdbqt \ --ligand lig_1.pdbqt \ --center_x 12.5 --center_y 8.3 --center_z 20.1 \ --size_x 22 --size_y 22 --size_z 22 \ --exhaustiveness 8 \ --num_modes 9 \ --energy_range 4 \ --cpu 4 \ --seed 42 \ --out lig_1_out.pdbqt \ --log lig_1.log--out是含构象的 PDBQT--log是文本日志里面第一张表就是各构象的亲和力。批量解析时读 log 比读 PDBQT 方便因为亲和力直接以文本形式列出来了。3.2 批量脚本遍历配体目录并调用 Vina下面这个脚本把配体目录里所有 PDBQT 逐个喂给 Vina参数从字典统一传入输出按配体名归档。import os import subprocess from pathlib import Path # 统一参数所有配体共用 RECEPTOR receptor.pdbqt CENTER (12.5, 8.3, 20.1) SIZE (22, 22, 22) EXHAUSTIVENESS 8 NUM_MODES 9 ENERGY_RANGE 4 CPU 4 SEED 42 ligand_dir Path(ligands_pdbqt) out_dir Path(docking_out) out_dir.mkdir(exist_okTrue) for lig in sorted(ligand_dir.glob(*.pdbqt)): name lig.stem out_pdbqt out_dir / f{name}_out.pdbqt log_file out_dir / f{name}.log cmd [ vina, --receptor, RECEPTOR, --ligand, str(lig), --center_x, str(CENTER[0]), --center_y, str(CENTER[1]), --center_z, str(CENTER[2]), --size_x, str(SIZE[0]), --size_y, str(SIZE[1]), --size_z, str(SIZE[2]), --exhaustiveness, str(EXHAUSTIVENESS), --num_modes, str(NUM_MODES), --energy_range, str(ENERGY_RANGE), --cpu, str(CPU), --seed, str(SEED), --out, str(out_pdbqt), --log, str(log_file), ] # 捕获异常单个失败不中断整批 try: subprocess.run(cmd, checkTrue, capture_outputTrue, textTrue) print(f[OK] {name}) except subprocess.CalledProcessError as e: print(f[FAIL] {name}: {e.stderr[:200]})subprocess.run的checkTrue让非零退出码抛异常配合 try 保证一个配体失败不会拖垮整批。capture_outputTrue把 Vina 的 stdout/stderr 收进变量失败时能打印前 200 字符定位问题。sorted()保证遍历顺序稳定方便对照结果。3.3 并行加速用 multiprocessing 把核吃满上面的串行脚本在 200 个配体、每个 30 秒的情况下要跑近两小时。机器有多核时用进程池并行能线性提速。注意 Vina 自己也有--cpu参数并行时要把单任务的--cpu调小避免核数超订。from multiprocessing import Pool import subprocess from pathlib import Path def run_one(lig_path): lig Path(lig_path) name lig.stem out_pdbqt fdocking_out/{name}_out.pdbqt log_file fdocking_out/{name}.log cmd [ vina, --receptor, receptor.pdbqt, --ligand, str(lig), --center_x, 12.5, --center_y, 8.3, --center_z, 20.1, --size_x, 22, --size_y, 22, --size_z, 22, --exhaustiveness, 8, --num_modes, 9, --cpu, 2, --seed, 42, --out, out_pdbqt, --log, log_file, ] try: subprocess.run(cmd, checkTrue, capture_outputTrue, textTrue) return (name, OK) except subprocess.CalledProcessError as e: return (name, fFAIL: {e.stderr[:100]}) if __name__ __main__: ligs [str(p) for p in sorted(Path(ligands_pdbqt).glob(*.pdbqt))] # 8 个进程并行每个 Vina 用 2 核共 16 核 with Pool(processes8) as pool: for name, status in pool.imap_unordered(run_one, ligs): print(name, status)Pool(processes8)开 8 个进程每个 Vina 用--cpu 2总占用 16 核。imap_unordered谁先跑完谁先返回比map更早拿到结果。进程数不是越多越好超过物理核数后上下文切换开销会吃掉收益一般设成物理核数 / 单任务cpu。注意并行写日志时如果多个进程写同一个文件会互相覆盖所以每个配体必须用独立的 log 文件名脚本里用name做前缀就是为了这个。4. 批量对接结果的解析、排序与筛选4.1 从 Vina log 里提取亲和力Vina 的 log 文件里有一张以-------------------------------------分隔的表第一列是 mode第二列是 affinitykcal/mol。用正则或按行切分都能提。import re from pathlib import Path def parse_vina_log(log_path): 返回 (最佳亲和力, 该构象的rmsd_lb, rmsd_ub) text Path(log_path).read_text() # 匹配表格数据行数字 空格 数字 空格 数字 空格 数字 pattern re.compile(r^\s*(\d)\s(-?\d\.\d)\s(\d\.\d)\s(\d\.\d), re.M) matches pattern.findall(text) if not matches: return None # 第一行就是最优构象 mode, affinity, rmsd_lb, rmsd_ub matches[0] return float(affinity), float(rmsd_lb), float(rmsd_ub) # 遍历所有 log汇总成列表 results [] for log in sorted(Path(docking_out).glob(*.log)): parsed parse_vina_log(log) if parsed: affinity, lb, ub parsed results.append((log.stem, affinity, lb, ub)) # 按亲和力升序越负结合越强 results.sort(keylambda x: x[1]) for name, aff, lb, ub in results[:10]: print(f{name}\t{aff}\t{lb}\t{ub})正则^\s*(\d)\s(-?\d\.\d)...抓的是表格数据行re.M让^匹配每行开头。matches[0]取第一行因为 Vina 按亲和力从优到劣排列第一行就是最佳构象。亲和力是负值越负表示结合越强所以排序用升序。4.2 结果汇总成 CSV 并做初筛把结果写成 CSV方便丢进 Excel 或 pandas 继续分析。初筛一般设一个亲和力阈值比如-7.0 kcal/mol再叠加配体效率亲和力除以重原子数做二次排序。import csv # 假设已有 results 列表和每个配体的重原子数字典 heavy_atoms with open(docking_summary.csv, w, newline) as f: writer csv.writer(f) writer.writerow([ligand, affinity, rmsd_lb, rmsd_ub, heavy_atoms, ligand_efficiency]) for name, aff, lb, ub in results: n heavy_atoms.get(name, 0) le aff / n if n else writer.writerow([name, aff, lb, ub, n, le])配体效率ligand efficiency是亲和力除以重原子数用来比较不同大小分子的结合质量。一个 30 重原子、亲和力 -9 的分子配体效率是 -0.3一个 15 重原子、亲和力 -6 的分子配体效率是 -0.4后者单位原子的结合贡献更高往往更值得优化。4.3 用 RMSD 判断构象是否收敛rmsd_lb和rmsd_ub是相对最优构象的下界和上界 RMSD。如果第二、第三构象的 RMSD 都很大说明搜索没收敛可能盒子设小了或--exhaustiveness不够。批量结果里如果大量配体的 top 构象 RMSD 分布很散就该回头检查盒子参数而不是继续往下筛。现象可能原因处理所有配体亲和力都接近 0盒子没盖住口袋重算盒子中心top 构象 RMSD 普遍偏大exhaustiveness 太低提到 16–32 重跑个别配体报错退出PDBQT 格式异常单独检查该配体转换结果和已知抑制剂对不上受体加氢/电荷有问题重做受体准备5. 批量对接的进阶技巧与可复现性保障5.1 固定随机种子让结果可复现Vina 的搜索带随机性同一个配体跑两次结果可能不同。--seed固定后相同输入和参数会给出相同输出。批量筛选里这一点很重要否则你没法解释为什么昨天排第一的配体今天掉到第五。我一般把 seed 写进配置所有配体共用同一个值保证整批结果内部一致。5.2 用配置文件替代长命令行参数一多命令行会变得又长又难维护。Vina 支持把参数写进配置文件用--config加载# vina_config.txt receptor receptor.pdbqt center_x 12.5 center_y 8.3 center_z 20.1 size_x 22 size_y 22 size_z 22 exhaustiveness 8 num_modes 9 energy_range 4 cpu 2 seed 42调用时只需vina --config vina_config.txt --ligand lig_1.pdbqt --out lig_1_out.pdbqt --log lig_1.log。批量脚本里配体、输出、日志三个路径是变量其余全从配置读脚本更干净改参数也只改一处。5.3 断点续跑跳过已完成的配体批量任务跑到一半中断是常事。与其从头再来不如在脚本开头检查输出文件是否已存在且非空存在就跳过。from pathlib import Path def already_done(name): log Path(fdocking_out/{name}.log) out Path(fdocking_out/{name}_out.pdbqt) # 两个文件都存在且 log 里有亲和力表格才算完成 if log.exists() and out.exists() and log.stat().st_size 0: return True return False # 在主循环里 for lig in sorted(ligand_dir.glob(*.pdbqt)): if already_done(lig.stem): print(f[SKIP] {lig.stem}) continue # ... 调用 Vina判断完成不能只看文件存在因为中断时可能生成了空文件。加上st_size 0能过滤掉大部分半成品。更严格的做法是检查 log 里是否含亲和力表格用前面那个正则匹配一下匹配到才算真完成。5.4 结果验证用已知配体做阳性对照批量流程跑通后别急着信结果。拿一个已知的活性配体比如共晶结构里的原配体重新对接一遍看 Vina 给出的最优构象和晶体构象的 RMSD 是否在 2 Å 以内。如果 RMSD 很大说明盒子、受体准备或参数有问题整批结果的可信度都要打问号。这个阳性对照是批量对接里最容易被跳过、也最不该跳过的一步。本文还有配套的精品资源点击获取 Autodock Vina批量分子对接：Python脚本实现与并行加速指南
简介这份资料面向生物信息学与计算化学方向的研究者、研究生及药物设计初学者聚焦于在Ubuntu 18.04环境下搭建Autodock Vina分子对接流程并结合Slurm调度器实现集群上的批量对接任务。内容涵盖依赖库安装、Open Babel与MGLTools配置、受体与配体预处理脚本的使用以及通过conf.txt控制CPU数量与搜索空间等关键参数最终输出pdbqt并转换为sdf便于后续分析同时详细说明Slurm的munge认证、slurmctld与slurmd服务启动及slurm.conf编写方法。资源包为1个PDF文件大小约1013KB以图文步骤形式呈现完整操作链路。目前已有5185人学习适合希望从单机对接进阶到集群批量计算、提升虚拟筛选效率的读者参考。1. 从一次 200 个配体的对接任务说起单个分子对接跑一次 Autodock Vina 只要几十秒很多人第一次接触时觉得这东西挺快。但真正落到虚拟筛选场景比如你手上有 200 个从天然产物库筛出来的候选配体要对着同一个靶点蛋白逐个跑手动改配置、改文件名、点运行、等结果一天就搭进去了。批量分子对接要解决的核心问题不是怎么对接而是怎么把重复劳动交给脚本同时保证每个配体的参数一致、结果可追溯。Autodock Vina 本身是命令行工具这个特性决定了它天然适合被批处理脚本调用。它接受--receptor、--ligand、--center_x/y/z、--size_x/y/z这些参数输出一个包含若干结合构象和亲和力的文本结果。批量化的关键就在于把配体列表、受体文件、对接盒子参数抽成配置用一个循环或并行框架去驱动 Vina再把散落的输出统一解析成表格。适合读这篇的人是已经会跑单次对接、现在要把规模从 1 提到 100 以上的生信或药物设计从业者。2. Autodock Vina 批量对接前的文件准备与参数确定2.1 受体与配体的格式转换从 PDB 到 PDBQTVina 不直接吃 PDB 文件受体和配体都得转成 PDBQT。受体转换用 MGLTools 里的prepare_receptor4.py配体用prepare_ligand4.py。这一步是批量流程里最容易翻车的地方因为配体来源五花八门有的带氢有的不带有的电荷没加。# 受体准备去水、加氢、加电荷输出 PDBQT prepare_receptor4.py -r receptor.pdb -o receptor.pdbqt \ -A hydrogens -U nphs_lps_waters # 单个配体准备 prepare_ligand4.py -l ligand.sdf -o ligand.pdbqt -A hydrogens-A hydrogens表示自动加氢-U nphs_lps_waters表示删除非极性氢、孤立水分子和孤对电子。批量场景下配体往往是 SDF 或 MOL2 格式需要先拆分再逐个转换。常见做法是用 Open Babel 先把多分子 SDF 拆成单文件# 把含多个分子的 SDF 拆成独立文件 obabel library.sdf -O lig_.mol2 -m-m是 split 模式会生成lig_1.mol2、lig_2.mol2这样的序列文件。拆完之后再写个循环批量转 PDBQT比一个个手点靠谱得多。2.2 对接盒子中心与尺寸怎么定才不白跑对接盒子grid box决定了搜索空间。盒子太小配体可能被截断在边缘构象不合理盒子太大搜索空间爆炸Vina 跑得慢还容易给出假阳性。定盒子的标准做法是围绕已知活性位点或共晶配体中心。如果你有共晶结构最省事的办法是算共晶配体的几何中心# 用 MDAnalysis 算共晶配体中心作为盒子中心 import MDAnalysis as mda u mda.Universe(complex.pdb) lig u.select_atoms(resname LIG) # 换成你的配体残基名 center lig.center_of_mass() print(center_x , round(center[0], 3)) print(center_y , round(center[1], 3)) print(center_z , round(center[2], 3))center_of_mass()返回质心坐标直接填进 Vina 的--center_x/y/z。盒子尺寸一般让每个方向比配体最大跨度多留 8–12 Å常见取值是 20×20×20 到 30×30×30。没有共晶结构时用 fpocket、SiteMap 之类的工具预测口袋或者干脆把整个蛋白包进去做盲对接但盲对接的假阳性率会明显上升。2.3 参数表批量对接必须统一的几个值批量对接最忌讳每个配体用不同参数结果没法横向比较。下面这张表是我一般会固定下来的参数写进配置文件里让所有配体共用。参数含义常用取值说明--exhaustiveness搜索彻底程度8–32批量筛选用 8精确重跑用 32--num_modes输出构象数9默认 9够看排名--energy_range构象能量差上限3–4单位 kcal/mol--cpu使用核数按机器定单任务别占满留给并行--seed随机种子固定值保证可复现--exhaustiveness是批量场景里最需要权衡的参数。它和运行时间近似线性8 和 32 的耗时差 4 倍左右但打分排名在前几位的构象往往差别不大。做初筛时用 8 快速过一遍挑出 top 10% 再用 32 精算是性价比比较高的策略。3. 用 Python 驱动 Vina 实现批量分子对接3.1 单次 Vina 命令的拆解先把单次命令写清楚批量才有模板可套。一条完整的 Vina 命令长这样vina --receptor receptor.pdbqt \ --ligand lig_1.pdbqt \ --center_x 12.5 --center_y 8.3 --center_z 20.1 \ --size_x 22 --size_y 22 --size_z 22 \ --exhaustiveness 8 \ --num_modes 9 \ --energy_range 4 \ --cpu 4 \ --seed 42 \ --out lig_1_out.pdbqt \ --log lig_1.log--out是含构象的 PDBQT--log是文本日志里面第一张表就是各构象的亲和力。批量解析时读 log 比读 PDBQT 方便因为亲和力直接以文本形式列出来了。3.2 批量脚本遍历配体目录并调用 Vina下面这个脚本把配体目录里所有 PDBQT 逐个喂给 Vina参数从字典统一传入输出按配体名归档。import os import subprocess from pathlib import Path # 统一参数所有配体共用 RECEPTOR receptor.pdbqt CENTER (12.5, 8.3, 20.1) SIZE (22, 22, 22) EXHAUSTIVENESS 8 NUM_MODES 9 ENERGY_RANGE 4 CPU 4 SEED 42 ligand_dir Path(ligands_pdbqt) out_dir Path(docking_out) out_dir.mkdir(exist_okTrue) for lig in sorted(ligand_dir.glob(*.pdbqt)): name lig.stem out_pdbqt out_dir / f{name}_out.pdbqt log_file out_dir / f{name}.log cmd [ vina, --receptor, RECEPTOR, --ligand, str(lig), --center_x, str(CENTER[0]), --center_y, str(CENTER[1]), --center_z, str(CENTER[2]), --size_x, str(SIZE[0]), --size_y, str(SIZE[1]), --size_z, str(SIZE[2]), --exhaustiveness, str(EXHAUSTIVENESS), --num_modes, str(NUM_MODES), --energy_range, str(ENERGY_RANGE), --cpu, str(CPU), --seed, str(SEED), --out, str(out_pdbqt), --log, str(log_file), ] # 捕获异常单个失败不中断整批 try: subprocess.run(cmd, checkTrue, capture_outputTrue, textTrue) print(f[OK] {name}) except subprocess.CalledProcessError as e: print(f[FAIL] {name}: {e.stderr[:200]})subprocess.run的checkTrue让非零退出码抛异常配合 try 保证一个配体失败不会拖垮整批。capture_outputTrue把 Vina 的 stdout/stderr 收进变量失败时能打印前 200 字符定位问题。sorted()保证遍历顺序稳定方便对照结果。3.3 并行加速用 multiprocessing 把核吃满上面的串行脚本在 200 个配体、每个 30 秒的情况下要跑近两小时。机器有多核时用进程池并行能线性提速。注意 Vina 自己也有--cpu参数并行时要把单任务的--cpu调小避免核数超订。from multiprocessing import Pool import subprocess from pathlib import Path def run_one(lig_path): lig Path(lig_path) name lig.stem out_pdbqt fdocking_out/{name}_out.pdbqt log_file fdocking_out/{name}.log cmd [ vina, --receptor, receptor.pdbqt, --ligand, str(lig), --center_x, 12.5, --center_y, 8.3, --center_z, 20.1, --size_x, 22, --size_y, 22, --size_z, 22, --exhaustiveness, 8, --num_modes, 9, --cpu, 2, --seed, 42, --out, out_pdbqt, --log, log_file, ] try: subprocess.run(cmd, checkTrue, capture_outputTrue, textTrue) return (name, OK) except subprocess.CalledProcessError as e: return (name, fFAIL: {e.stderr[:100]}) if __name__ __main__: ligs [str(p) for p in sorted(Path(ligands_pdbqt).glob(*.pdbqt))] # 8 个进程并行每个 Vina 用 2 核共 16 核 with Pool(processes8) as pool: for name, status in pool.imap_unordered(run_one, ligs): print(name, status)Pool(processes8)开 8 个进程每个 Vina 用--cpu 2总占用 16 核。imap_unordered谁先跑完谁先返回比map更早拿到结果。进程数不是越多越好超过物理核数后上下文切换开销会吃掉收益一般设成物理核数 / 单任务cpu。注意并行写日志时如果多个进程写同一个文件会互相覆盖所以每个配体必须用独立的 log 文件名脚本里用name做前缀就是为了这个。4. 批量对接结果的解析、排序与筛选4.1 从 Vina log 里提取亲和力Vina 的 log 文件里有一张以-------------------------------------分隔的表第一列是 mode第二列是 affinitykcal/mol。用正则或按行切分都能提。import re from pathlib import Path def parse_vina_log(log_path): 返回 (最佳亲和力, 该构象的rmsd_lb, rmsd_ub) text Path(log_path).read_text() # 匹配表格数据行数字 空格 数字 空格 数字 空格 数字 pattern re.compile(r^\s*(\d)\s(-?\d\.\d)\s(\d\.\d)\s(\d\.\d), re.M) matches pattern.findall(text) if not matches: return None # 第一行就是最优构象 mode, affinity, rmsd_lb, rmsd_ub matches[0] return float(affinity), float(rmsd_lb), float(rmsd_ub) # 遍历所有 log汇总成列表 results [] for log in sorted(Path(docking_out).glob(*.log)): parsed parse_vina_log(log) if parsed: affinity, lb, ub parsed results.append((log.stem, affinity, lb, ub)) # 按亲和力升序越负结合越强 results.sort(keylambda x: x[1]) for name, aff, lb, ub in results[:10]: print(f{name}\t{aff}\t{lb}\t{ub})正则^\s*(\d)\s(-?\d\.\d)...抓的是表格数据行re.M让^匹配每行开头。matches[0]取第一行因为 Vina 按亲和力从优到劣排列第一行就是最佳构象。亲和力是负值越负表示结合越强所以排序用升序。4.2 结果汇总成 CSV 并做初筛把结果写成 CSV方便丢进 Excel 或 pandas 继续分析。初筛一般设一个亲和力阈值比如-7.0 kcal/mol再叠加配体效率亲和力除以重原子数做二次排序。import csv # 假设已有 results 列表和每个配体的重原子数字典 heavy_atoms with open(docking_summary.csv, w, newline) as f: writer csv.writer(f) writer.writerow([ligand, affinity, rmsd_lb, rmsd_ub, heavy_atoms, ligand_efficiency]) for name, aff, lb, ub in results: n heavy_atoms.get(name, 0) le aff / n if n else writer.writerow([name, aff, lb, ub, n, le])配体效率ligand efficiency是亲和力除以重原子数用来比较不同大小分子的结合质量。一个 30 重原子、亲和力 -9 的分子配体效率是 -0.3一个 15 重原子、亲和力 -6 的分子配体效率是 -0.4后者单位原子的结合贡献更高往往更值得优化。4.3 用 RMSD 判断构象是否收敛rmsd_lb和rmsd_ub是相对最优构象的下界和上界 RMSD。如果第二、第三构象的 RMSD 都很大说明搜索没收敛可能盒子设小了或--exhaustiveness不够。批量结果里如果大量配体的 top 构象 RMSD 分布很散就该回头检查盒子参数而不是继续往下筛。现象可能原因处理所有配体亲和力都接近 0盒子没盖住口袋重算盒子中心top 构象 RMSD 普遍偏大exhaustiveness 太低提到 16–32 重跑个别配体报错退出PDBQT 格式异常单独检查该配体转换结果和已知抑制剂对不上受体加氢/电荷有问题重做受体准备5. 批量对接的进阶技巧与可复现性保障5.1 固定随机种子让结果可复现Vina 的搜索带随机性同一个配体跑两次结果可能不同。--seed固定后相同输入和参数会给出相同输出。批量筛选里这一点很重要否则你没法解释为什么昨天排第一的配体今天掉到第五。我一般把 seed 写进配置所有配体共用同一个值保证整批结果内部一致。5.2 用配置文件替代长命令行参数一多命令行会变得又长又难维护。Vina 支持把参数写进配置文件用--config加载# vina_config.txt receptor receptor.pdbqt center_x 12.5 center_y 8.3 center_z 20.1 size_x 22 size_y 22 size_z 22 exhaustiveness 8 num_modes 9 energy_range 4 cpu 2 seed 42调用时只需vina --config vina_config.txt --ligand lig_1.pdbqt --out lig_1_out.pdbqt --log lig_1.log。批量脚本里配体、输出、日志三个路径是变量其余全从配置读脚本更干净改参数也只改一处。5.3 断点续跑跳过已完成的配体批量任务跑到一半中断是常事。与其从头再来不如在脚本开头检查输出文件是否已存在且非空存在就跳过。from pathlib import Path def already_done(name): log Path(fdocking_out/{name}.log) out Path(fdocking_out/{name}_out.pdbqt) # 两个文件都存在且 log 里有亲和力表格才算完成 if log.exists() and out.exists() and log.stat().st_size 0: return True return False # 在主循环里 for lig in sorted(ligand_dir.glob(*.pdbqt)): if already_done(lig.stem): print(f[SKIP] {lig.stem}) continue # ... 调用 Vina判断完成不能只看文件存在因为中断时可能生成了空文件。加上st_size 0能过滤掉大部分半成品。更严格的做法是检查 log 里是否含亲和力表格用前面那个正则匹配一下匹配到才算真完成。5.4 结果验证用已知配体做阳性对照批量流程跑通后别急着信结果。拿一个已知的活性配体比如共晶结构里的原配体重新对接一遍看 Vina 给出的最优构象和晶体构象的 RMSD 是否在 2 Å 以内。如果 RMSD 很大说明盒子、受体准备或参数有问题整批结果的可信度都要打问号。这个阳性对照是批量对接里最容易被跳过、也最不该跳过的一步。本文还有配套的精品资源点击获取