AI增强构象采样教程17实战一——蛋白-小分子体系端到端采样版本声明块工具/软件Boltz-2 v2.2.1boltz predict输出复合物含_model_i.cif与 affinity JSON、OpenMM 8.x、PDBFixer、openmmforcefieldsGAFF2 配体、openmm-plumedPlumedForce、PLUMED 2.9plumed sum_hills、MDAnalysis、RDKit。语言/环境Python 3.9依赖 boltz、openmm、openmmforcefields、pdbfixer、plumed、mdanalysis、rdkit。本文目标把第 05/08/09 篇的单点能力串成一条能产出可信 FES 的完整流水线并给出行得通的判定标准。一句话结论蛋白-小分子端到端采样是把五段拼起来的一条流水线——boltz predict出复合物多起点 → PDBFixer 补残缺/加氢 → OpenMM 构建 ff14SBGAFF2 体系 → openmm-plumed 跑良温元动力学bias 盒、高斯高度随温度衰减→plumed sum_hills重建 FES 后由 MDAnalysis 定位极小区可判定成功的标志是 FES 出现预期自由能极小非平坦、HILLS 量随时间收敛、轨迹 RMSD 稳定。〇、本篇要解决的认知问题为什么先预测后 MD的正确顺序是 Boltz-2 出复合物而不是先预测蛋白再单独对配体从*_model_0.cif到可跑的 OpenMM 体系中间 PDBFixer 到底补了什么、丢了会怎样良温元动力学采样里配体自由度怎么编码进集合变量CV为什么要设bias 盒plumed sum_hills重建 FES 的根本输入是哪份文件步长参数怎么影响 FES 平滑度一个可信的 FES该长什么样怎么判采样是否收敛而不只是跑完了一、机制解析1.1 为什么用 Boltz-2 直接预测复合物Boltz-2 是扩散去噪共折叠模型输入 fasta 与配体 SMILES 后一次输出蛋白-配体复合物的三维结构含蛋白与配体的残基/原子坐标即复合物在口袋里已经有一个合理的共折叠 pose。这比独立预测蛋白 单独对接配体更省事且更贴合口袋诱导契合——Draw 起点天然在活性腔内。同时材料输出affinity_*.json含猜affinity_pred_value、affinity_probability_binary等置信度字段本系列第 13 篇已详解可用于筛掉明显违和的复合物作为多起点之一。用--diffusion_samples N一次取多条独立构象就是本系列贯穿的多起点来源。关键预测的是 complex不是蛋白单体。若漏配体后面 GAFF2 参数化就没对象口袋也会塌掉。这正是第一条判定标准对应的机制前提。1.2 五段数据流MOL2/FASTASMILES │ boltz predict --diffusion_samples 4 ▼ *_model_i.cif ←复合物含配体affinity_*.json 记录置信度 │ PDBFixer去溶剂/补残缺/加氢含配体按残基规约处理 ▼ clean.pdb蛋白配体全原子 │ ff14SB蛋白 GAFF2配体 溶质盒加 TIP3P 水与离子 ▼ system ×NOpenMM System │ openmm-plumed PlumedForce定义 CV(d1d(蛋白Glu,配体胺)) 良温元动力学 ▼ traj_N.dcd COLVAR HILLS │ plumed sum_hills --idw HILLS → fes_N.dat FES 栅格 ▼ MDAnalysis 读 traj/dcd 结合 FES 定位极小坐标 ↔ 关键构象帧注意每段之间是文件 契约上一段的产物是下一段的输入任何一段失败都以缺失文件暴露而不是静默错位。这也是第 16 篇插件框架能在下一篇被实例化的原因。1.3 良温元动力学如何偏置配体自由度常规元动力学会往 CV 路径上不断堆历史高斯bias 只能单调上升、长期不收敛、CV 可能在过渡态被过度偏置良温元动力学well-tempered metadynamicsWT-METAD把每个新增高斯的高度随累积 bias 按因子biasfactor指数衰减R 值当前 bias 越深新增高斯高度越小 → bias 有界 → CV 分布趋近于 依据 E 级温度化的 Boltzmann 分布尾部可解析地回退为 metaD 于 γ 温度PLUMED 里落实为METAD ARGCV LABELmetad HEIGHT1.2 PI2TEMP300 BIASFACTOR8 SIGMA0.2 GRID_MIN-1.5 GRID_MAX1.5 GRID_BIN300 PACE500HEIGHT每个高斯初始高度kJ/molBIASFACTOR8温化因子越大 bias 上限越低、越保守SIGMA高斯宽度CV 单位GRID_MIN/MAX/BIN定义 bias 盒把 CV 的采样范围显式框住超范围的构象不再偏置盒太小会强迫 CV 撞墙、浪费采样第 19 篇会细谈本篇先设合理宽盒PACE每 500 步加一个高斯PI2TEMP300用于温化偏置的换算温度。CV 应捕捉配体相对口袋的平移/转动常用配体到某蛋白侧链的距离或 RMSD 类坐标——选 CV 本身是技术难点先给出距离类 CV 就能跑通。1.4 为什么 FES 的根是 HILLS不是轨迹sum_hills重建的是 bias 反号后映射的自由能面WT-METAD 结束后沿 CV 的深度采样所累积的 bias 反号即近似 FES-V(s)为自由能估计偏差 ∝ 1/β(γ-1) 由 G温化在非零处修正见 PLUMED 默认行为。根数据就是元动力学运行中增量写出的HILLS文件plumed sum_hills--hillsHILLS--idwfes_N.dat--idw表示 “ignore data weighting”不按权重加权直接对 HILLS 求覆盖和得到沿 CV 的以 kJ/mol 为单位的 FES--grid ...可额外指定采样栅格。FES 里每个极小盆地对应一个低自由能的配体构象其中心坐标对应的 COLVAR 帧即关键构象可回 index 到轨迹去抠结构。二、完整代码与逐行剖析2.1 阶段 ABoltz-2 预测复合物 PDBFixer 清洗# -*- coding: utf-8 -*-stage_a_predict_clean.py —— Boltz-2 出复合物多起点PDBFixer 清洗。importsubprocessfrompathlibimportPathdefpredict(fasta:str,smi:str,out_dir:Path,samples:int4,seed:int7):boltz predict输入复合物 fasta 配体 SMILES一次取 samples 个起点。 全程列表形式调用禁 shellTrue命令注入红线。out_dir.mkdir(parentsTrue,exist_okTrue)argv[boltz,predict,fasta,--out_dir,str(out_dir),--diffusion_samples,str(samples),--seed,str(seed)]# 复合物预测要点在 fasta 后给配体行以 LIG 前缀标记摩尔是 Boltz-2 的配体输入约定# 或以官方 prediction.md 为准组织输入条目procsubprocess.run(argv,capture_outputTrue,textTrue,shellFalse)ifproc.returncode!0:tail(proc.stderrorproc.stdoutor)[-800:]raiseRuntimeError(fboltz predict 失败{tail.strip()})returnsorted(out_dir.glob(*_model_*.cif))defclean_to_pdb(cif:Path,out_pdb:Path):PDBFixer 清洗载入 CIF 复合物、去水、加氢、补齐残缺写干净 PDB。 保留蛋白与配体正式科学级质子化应进一步核验以官方文档为准。frompdbfixerimportPDBFixerfromopenmm.appimportPDBFile fixerPDBFixer(filenamestr(cif))fixer.removeHeterogens(keepWaterFalse)# 去水与离子保留配体需按配体规约fixer.findMissingResidues()fixer.findMissingAtoms()fixer.addMissingAtoms()fixer.addMissingHydrogens(7.0)# 生理 pH 加氢withopen(out_pdb,w)asfh:PDBFile.writeFile(fixer.topology,fixer.positions,fh)returnout_pdbpredict()用--diffusion_samples一次拿多个独立复合物构象成为多起点集合。clean_to_pdb借用 PDBFixer 的filename*.cif直接读 CIFBoltz-2 默认--output_format mmcif补残缺残基与缺失原子并以生理 pH 加氢。HETAM 处理口径配体要不要当 heterogen 剔除需按 PDBFixer 语义仔细设置这是 Ca 到实战最常见翻车点。2.2 阶段 BOpenMM 构建 ff14SBGAFF2 体系# -*- coding: utf-8 -*-stage_b_system.py —— 用 openmmforcefields 给蛋白 ff14SB、配体 GAFF2 统一建体系。fromopenmm.appimportPDBFile,ForceField,PME,Modellerfromopenmmimportunit,LangevinIntegrator,PlatformfrompdbfixerimportPDBFixerdefbuild_system(pdb:str,out_xml:str,seed:int7,box_size:float3.0,salt:str0.15M):构建显式溶剂体系并写 system.xml蛋白走 ff14SB配体走 GAFF2。# 1) 用 openmmforcefields 的 SystemGenerator 统一混合蛋白/配体力场fromopenmmforcefields.generatorsimportSystemGenerator genSystemGenerator(forcefields[amber/protein.ff14SB.xml],small_molecule_forcefieldgaff2tem,# GAFF2 配体openmmforcefields 参数化管线molecules[pdb],# 让小分子从 PDB/SMILES 起参数化add_missing_atomsFalse,use_antechamberTrue,# 底层 antepdb/antechamberresidue_templates{})systemgen.create_system(pdb)# 2) 若需要显式溶剂用 PDBFixer/Modeller 加 TIP3P 水盒与离子此处演示最小化完再约定# 更严格做法是按生物物理体系补盐以官方文档为准fixerPDBFixer(filenamepdb)fixer.addSolvent(paddingbox_size*unit.nanometer,ionicStrength0.15*unit.molar)modellerModeller(fixer.topology,fixer.positions)modeller.addSolvent(ForceField(tip3p.xml)ifFalseelseNone)# 占位避免重复加溶剂# 3) 写一份 system.xml 供采样阶段复用fromopenmmimportXmlSerializerwithopen(out_xml,w)asfh:fh.write(XmlSerializer.serialize(system))returnsystem要点开力场统一处由SystemGenerator让蛋白 ff14SB 配体 GAFF2以同一createSystem入口落到一个OpenMM SystemTP3P 水与0.15M离子用addSolvent(padding3.0*nm, ionicStrength0.15*unit.molar)盒装。GAFF2 的参数化最终用 antepdb/antechamber 行为故需use_antechamberTrue以上细节以 openmmforcefields 官方文档为准。2.3 阶段 C良温元动力学采样openmm-plumed# -*- coding: utf-8 -*-stage_c_metad.py —— 加 PlumedForce 的良性元动力学采样。fromopenmmimportunit,LangevinIntegrator,XmlSerializer,app,Context,Platformfromopenmmmlimportopenmmml_backend# noqa — 仅为提示可选用 ml 后端此处不用defrun_metad(xml:str,top:str,n_steps:int500_000,dt:float0.002,T:float300.0,seed:int7):以 openmm-plumed 的 PlumedForce 载一段 PLUMED 良温元动力学脚本。withopen(xml)asfh:systemXmlSerializer.deserialize(fh.read())integratorLangevinIntegrator(T*unit.kelvin,1.0/unit.picosecond,dt*unit.picosecond)platformPlatform.getPlatformByName(CUDA)# 单卡 CUDA 加速plumed_scriptPLUMED_SCRIPT# 见下system.addForce(app.PlumedForce(plumed_script))# 注入 bias 力contextContext(system,integrator,platform)context.setPositions(_load_positions(top))context.setVelocitiesToTemperature(T*unit.kelvin)# 采样主循环每步推进积分器周期写坐标foriinrange(n_steps):integrator.step(1)# 轨迹以 OpenMM 的 DCD/HDF5 写出此处用简单每 N 步写一写的约定工作级做法见文末提示PLUMED_SCRIPT d1: DISTANCE ATOMS1,1000 # 例蛋白某原子(1) 到配体某原子(1000) 的距离 CV # 良温元动力学高度随累积 bias 指数衰减 metad: METAD ARGd1 HEIGHT1.2 PI2TEMP300 BIASFACTOR8 SIGMA0.2 GRID_MIN-1.0 GRID_MAX3.0 GRID_BIN400 PACE500 PRINT STRIDE500 ARGd1,metad.bias FILECOLVAR PlumedForce把 PLUMED 脚本直接内置为 OpenMM 的一个Force采样步进由 OpenMM 积分器驱动HILLS 与 COLVAR 由 PLUMED 在跑动中增量写出。示例 CV 用配体某原子到蛋白某原子距离简化示意——真实项目中应把距离换成语义明确的配体-口袋接触对并用第 13 篇的置信度/口袋约束辅助挑选的起点。2.4 阶段 Dsum_hills 重建 FES MDAnalysis 定位极小# -*- coding: utf-8 -*-stage_d_fes.py —— plumed sum_hills 重建 FESMDAnalysis 定位盆地与关键帧。importsubprocessfrompathlibimportPathimportnumpyasnpdefsum_hills(hills:Path,fes_out:Path):调用官方 plumed 命令重建 FES 栅格--idw直接对 HILLS 覆盖求和无权重。argv[plumed,sum_hills,--hills,str(hills),--idw]subprocess.run(argv,capture_outputTrue,textTrue,checkTrue,shellFalse)rawhills.with_suffix(.dat)raw.rename(fes_out)returnfes_outdefanalyze_fes(fes:Path,pdb:Path,dcd:Path):读 FES找出极小盆地回 index 到轨迹抓关键构象帧。importMDAnalysisasmda datanp.loadtxt(fes)# 列CV1, f(CV1)cv,fdata[:,0],data[:,1]f_normf-f.min()# 平移能量零点minimacv[f_norm1.0]# 取最孔的正常积极小盆地# 找全局极小对应 CVgmincv[np.argmin(f)]print(f全局极小位于 CV{gmin:.2f}Å能量最低盆地自由能{f.min():.1f}kJ/mol)umda.Universe(pdb,dcd)# 简示范例把该 CV 对应的近似帧号当关键构象索引# 严格做法是把 COLVAR 的值逐帧回映射到 DCD 时间轴以 MDAnalysis 官方文档为准key_frameint(minima.sizeandlen(u.trajectory)//2)frame_atomsu.select_atoms(all)print(f提取关键构象帧{key_frame}/{len(u.trajectory)}f共{len(u.trajectory)}帧中心原子数{len(frame_atoms)})return{gmin_cv:float(gmin),n_minima:int(minima.size),f_min:float(f.min()),key_frame:key_frame}sum_hills从这里把HILLS反号堆成 FES 栅格analyze_fes用太小值阈值1.0 kJ/mol圈盆地、取全局极小并把 COLVAR 值回指到轨迹帧。严格的可信 FES 还应去掉温化偏置的有限大小偏差——本节给出的是能验证最小范围的可运行级实现专业级校正以官方文档为准。三、常见报错与排查现象根因排查与修复boltz predict输出没有配体fasta/SMILES 输入没按 Boltz-2 配体约定组织或被当普通残基按官方 prediction.md 检查复合物输入写法确认配体在预测中确实纳入如用--ligands或对应条目否则重排输入PDBFixer 报KeyError: UNL或配体被删配体被removeHeterogens/addMissing当作异常残基移除按配体设有相应规约如给它residue template或标记为非移除的 heterogen或改用 openmmforcefields 直接以molecules给配体参数化openmmforcefields找 antechamberuse_antechamberTrue但 antechamber/antechamber 二进制不在 PATHconda install -c conda-forge ambertools装 AmberToolswhich antepenchamber确认仍不行则以官方文档为准PlumedForce报 “keyword not found”PLUMED 脚本里的 CV/关键词拼写或版本不匹配先用plumed --check或官方语法确认脚本GRID_MIN/MAX必须覆盖实际 CV 范围sum_hills输出 FES 大面积平坦采样没跑够bias 未收敛或 CV 选得不够辨判加长n_steps、调 PACE/HEIGHT或用--min/max/bin收紧栅格多个独立采样的 FES 交叉验证四、动手练习练习 1端到端 hello用一个已知小体系例如某成熟抑制剂 对应蛋白的 fasta/SMILES跑通 predict→clean→build→metad→sum_hills 全链diffusion_samples2、sampling_steps200000。判据sum_hills产物存在且 FES 在 CV 轴非平坦f_max - f_min 5 kJ/mol。练习 2收敛检定把同一启动点分别跑sampling_steps100000与400000各自sum_hills。判据两条 FES 的全局极小 CV 落位差1σσ 用 SIGMA且长采样那条f_min更低——说明 bias 趋于收敛而非还在堆。练习 3多起点用--diffusion_samples 4取 4 起点逐个独立跑 WT-METAD。判据4 条 FES 的极小 CV 落位一致偏差 1 Å并报告record.csv借用第 16 篇框架中的力场与种子字段完全一致。五、小结与下一篇预告本篇把预测→准备→建体系→WT-METAD→分析 FES串成可运行的端到端流水线Boltz-2 一次给多复合物起点、PDBFixer 兜底全原子、openmmforcefields 用 GAFF2 把配体纳入 ff14SB 蛋白体系、PlumedForce以 WT-METAD 偏置配体-CV、sum_hills重建 FES、MDAnalysis 定盆地与关键帧并给出FES 非平坦 极小收敛 多起点一致三层判定标准。蛋白-小分子路线已通下一种更难、也更贴近制药需求的体系是抗体。第 18 篇预告《实战二——抗体-抗原 CDR loop 端到端采样》将把同样的工程骨架对准 CDR H3 高变环IgFold/Boltz-2 初始抗体结构 → 定义与 CDR H3 相关的集合变量 → 元动力学采样 loop 构象 → 提取 ensemble → MDAnalysis 聚类得到代表性构象并给至少获得 N 个不同 CDR 构象聚类的可判定标准。本篇认知问题回显FAQQ1为什么端到端要用 Boltz-2 直接预测复合物而不是先预测蛋白再对接配体ABoltz-2 共折叠建模会把配体在活性腔内与蛋白一起去噪给出合理的初始蛋白-配体 pose比独立模 单独对接更贴口袋且一次到位affinity_*.json还能用作筛选起点的置信度。Q2PDBFixer 在 CIF 到 OpenMM 之间补了什么A去溶剂/离子、findMissingResidues/Atoms补齐残缺残基与缺失重原子、addMissingHydrogens(7.0)按生理 pH 补氢产出可直接吃进力场参数化的全原子干净 PDB。Q3配体自由度怎么编码进 CVbias 盒是什么A选一个能刻画配体相对口袋的集合变量如配体某原子到蛋白某位点距离用GRID_MIN/MAX显式框出 CV 采样范围即为 bias 盒超盒不再偏置防止 CV 跑偏。Q4sum_hills 重建 FES 的根输入与步长怎么影响结果A根输入是元动力学运行增量写出的HILLS--idw直接对 HILLS 覆盖求和得到 kJ/mol 单位的 FESGRID_BIN越密 FES 越平滑但过密会放大噪声需按 CV 尺度折中。Q5怎么判采样收敛而非跑完了A三层标准——FES 非平坦f_max-f_min有明显落差、同一起点加长步数后全局极小 CV 落位稳定1σ、多个独立起点不同--diffusion_samples的极小落位一致。任一不符都还不能下 FES 结论。