4、受体文件准备:从PDB到PDBQT的转换流程(MGLTools / OpenBabel)
这是个老生常谈,却总有人翻车的话题。我刚开始做对接那会儿,也天真地以为,直接从PDB数据库拉个结构就能扔给Vina跑。结果呢?跑出来的打分我自己都不信。说白了,Vina不认原始的PDB格式。它要的是PDBQT——一个包含了原子类型、电荷和可旋转键信息的“智能”版本。
这一步没做好,你后面的所有努力都白费。复现性?不存在的。今天我们就把这个流程拆开揉碎了讲。
4.1 为什么非转不可?PDB文件缺什么
PDB文件记录了原子的三维坐标,但这对Vina来说远远不够。Vina需要知道这个原子是啥“类型”——不是元素符号,而是力场里的角色。比如碳原子,在Vina里分好几种:脂肪碳(C)、芳香碳(A)、碳氢化合物里的碳等等。PDB文件里可没这些信息。
电荷也是个问题。原始PDB通常不带电荷,或者电荷分布不完整。Vina的评分函数依赖静电势,没有电荷,计算就是瞎猜。我个人习惯,在转换前先用PyMOL看一眼结构,晶体水、配体、金属离子有没有残留。这些杂质会影响电荷计算,甚至让Vina崩溃。
4.2 两条路:MGLTools的ADT vs. OpenBabel
走哪条路,取决于你的工作流。我一般在Windows上搞可视化就开ADT,在Linux服务器上跑批量就用OpenBabel。两条路的结果基本一致,但细节上有些差异。
4.2.1 使用MGLTools ADT(带图形界面)
ADT(AutoDock Tools)是MGLTools套件里的一个组件。它把PDB到PDBQT的转换做成了一个向导式流程。
基本操作是这样:
- 启动ADT,点击 Grid → Macromolecule → Open,导入你的受体PDB文件。
- ADT自动弹出提示,询问是否“合并非极性氢”和“分配Gasteiger电荷”。都选Yes。
- 程序自动运行,生成
receptor.pdbqt,保存在当前目录。
嗯,这里要注意。ADT的默认设置是把水分子和配体都删掉,只保留蛋白质残基和金属离子。如果你需要保留水分子,必须在导入前手动编辑PDB文件。我在项目中遇到过这种需求——某些水分子是催化活性所必需的。那不保留水分子,结果就丢了。我的做法是在PDB文件里只保留我需要的水分子(通常是离活性中心3Å以内的),然后强行指定它们的“原子类型”为W。
4.2.2 使用OpenBabel(命令行)
OpenBabel是个化学文件格式转换瑞士军刀。它支持的格式比ADT多得多,速度也快得多。命令行一转,一秒搞定。
基本命令长这样:
obabel receptor.pdb -O receptor.pdbqt -h --partialcharge gasteiger
解释一下参数:
-h:为受体添加氢原子(通常加pH 7.4条件下的氢,中性)--partialcharge gasteiger:分配Gasteiger部分电荷。这是AutoDock力场默认用的电荷模型。
如果你有多个PDB文件要批量处理,加个循环:
for pdb in *.pdb; do
obabel "$pdb" -O "${pdb%.pdb}.pdbqt" -h --partialcharge gasteiger
done
4.3 PDBQT文件里到底有什么?
为了不让你觉得自己在盲转,我们拆开一个PDBQT看看。打开PDBQT文件,每行看起来像这样:
ATOM 1 N VAL A 1 12.345 23.456 34.567 1.00 0.00 0.000 N
ATOM 2 CA VAL A 1 13.456 24.567 35.678 1.00 0.00 0.000 C
ATOM 3 C VAL A 1 14.567 25.678 36.789 1.00 0.00 0.000 C
ATOM 4 O VAL A 1 15.678 26.789 37.890 1.00 0.00 0.000 OA
...
重点看最后一列:原子类型(AutoDock Type)。这决定了Vina怎么识别这个原子、怎么计分。常见的AutoDock原子类型有:
| 类型代码 | 描述 | 典型原子 |
|---|---|---|
| C | 脂肪碳(sp3碳) | CH3, CH2, CH, C |
| A | 芳香碳(sp2碳) | 苯环上的C |
| N | 氮原子(胺/酰胺) | N, NH, NH2 |
| OA | 氧原子(羰基/羟基) | C=O, O-H |
| HD | 极性氢(与N/O/S相连的H) | OH, NH中的H |
| NA | 氮原子(芳香环中的N) | 吡啶上的N |
除了原子类型,每个原子的部分电荷(Partial Charge)也是关键。在ATOM记录里,它是第9列(0.000那个位置)。Vina用它来计算静电相互作用。电荷分配不均,静电势就乱了。这就是为什么我反复强调Gasteiger电荷——这是AutoDock力场设计之初就配套的电荷模型,它能保证计算结果与实验数据的偏差最小。
4.4 自动加氢的陷阱
不管是ADT还是OpenBabel,都默认帮你“加氢”。但加氢可不是简单的给每个重原子补齐H原子数。它要考虑原子周围的化学环境、pH、甚至配位状态。你想想看,一个组氨酸(His)在生理pH下,它的咪唑环上到底哪个N是质子化的?加错了,电荷分布就全错了。
我踩过坑。有一回我处理一个含Zn²⁺离子的金属酶,用OpenBabel默认加氢,结果把Zn²⁺周围的配位水分子也加了H,变成了H₃O⁺。Vina跑的时候直接报计算错误。后来我手动删了那些水,重新加氢,问题解决。
4.5 非极性氢合并:为什么和多说一点点
PDBQT里的一个“原子”不一定对应一个真实原子。AutoDock有一套“非极性氢合并”机制。说白了,就是把连在碳上的氢(非极性氢)合并到它们连的碳原子上。这样做的目的是减少原子数量,提高计算效率。合并后,碳原子会带一个“伪原子”,它的电荷和范德华参数会综合考虑原本碳和它上面所有氢原子的贡献。
这个操作是必须的。在PDBQT文件的ROOT块里,你会看到类似这样的记录:
ROOT
ATOM 1 C UNK 1 1.234 2.345 3.456 1.00 0.00 0.000 C
ATOM 2 C UNK 1 1.987 3.210 4.321 1.00 0.00 0.000 C
ATOM 3 C UNK 1 1.111 2.222 3.333 1.00 0.00 0.000 C
ENDROOT
BRANCH 1 2
ATOM 4 C UNK 1 1.234 2.345 3.456 1.00 0.00 0.000 C
ENDBRANCH
注意看,ROOT部分都是“原”原子,但数量少了——碳原子的坐标和电荷都调整过了。这就是合并后的效果。
4.6 金属离子的处理
如果你受体里有金属离子(Zn²⁺、Mg²⁺、Fe²⁺等),处理要特别小心。AutoDock力场里没有定义金属离子的原子类型。常用的做法是:
- 保留金属离子,但指定它的原子类型为“金属”。在PDBQT里,金属离子的原子类型通常设为Mn(M代表金属,n是原子序数)。比如Zn就是Zn。
- 或者,保留金属离子但不给它电荷,把它当成一个“硬球”。这种方法会丢失静电信息,但至少Vina不会报错。
我记得有一次做个含Fe²⁺的氧化酶,我选择了第一种方法——给Fe²⁺指定Fe类型,并分配+2电荷。但Vina不认这个Fe类型,直接跳过所有相互作用。结果对接打分完全不靠谱。后来我查阅文献,发现AutoDock力场建议对过渡金属使用Lennard-Jones参数拟合。说白了,得手动在PDBQT里改参数,或者干脆把金属离子当成“刚性分子”的一部分来处理。
4.7 最后的一步检查
转完PDBQT,别急着跑Vina。养成检查文件的好习惯。我一般用这个命令看PDBQT的文件头:
head -20 receptor.pdbqt
至少要确保:
- 文件里有
ROOT和ENDROOT块。 - 原子数、残基数和你的PDB一致(加氢后会多几十个原子,但重原子数不变)。
- 每个原子都有非零的部分电荷(Gasteiger电荷通常不会是零)。
- 没有水分子残留(除非你刻意保留)。
检查完,再继续下一步。