4. PDBQT 文件格式详解

好，我们正式开始拆解 PDBQT 格式。说实话，我第一次接触这个格式的时候，觉得它跟 PDB 长得挺像，就没太在意。结果呢？在对接计算里栽了个跟头——ATOM 记录里的电荷和类型没写对，整个打分函数全跑偏了。从那以后，我养成了一个习惯：但凡涉及格式转换，先把目标格式的每一行、每个字段啃透。

PDBQT 说白了就是 PDB 的「增强版」。AutoDock 系列软件（包括 AutoDock Vina）不认识标准的 PDB，它只认 PDBQT。你可以把 PDBQT 理解成一套「定制化卡片」，它除了记录原子坐标，还额外塞了原子类型、部分电荷、以及配体内部的可旋转键信息。没有这些，对接程序根本不知道哪个键能转、哪个原子带多少电。

4.1 记录结构一览

PDBQT 文件的行记录主要分几大类。我建议你把下表当作「速查手册」，写脚本的时候贴在旁边。

记录类型	作用	常见示例
ATOM	记录原子坐标、类型、电荷	ATOM 1 C UNL 1 1.234 5.678 9.012 0.00 0.00 C 0.1234
HETATM	非标准残基原子（配体常用）	HETATM 2 O UNL 1 4.321 6.789 1.098 0.00 0.00 OA -0.5678
ROOT	标记配体的刚性根部分	ROOT 开头，ENDROOT 结束
BRANCH	标记可旋转分支的起点	BRANCH 3 6 表示原子3和6之间有关节
END	文件结束	END

嗯，这里要特别注意：如果是配体文件，通常用 HETATM 开头，而不是 ATOM。受体蛋白文件才用 ATOM。我做格式检测的时候，第一件事就是 grep 一下文件头，看看用的是哪种记录。

4.2 核心字段深度解析

4.2.1 原子名称与坐标

这部分跟 PDB 基本一致：ATOM 或 HETATM 开头，后面跟着原子序号、原子名、残基名、链 ID、残基号、以及 X/Y/Z 坐标。注意，PDBQT 默认不写 B 因子和占有率，或者直接写 0.00。

HETATM 1 C3 UNL 1  -0.159  -0.044   1.059  0.00  0.00   C  -0.1325

你看最后两位：C 是原子类型，-0.1325 是部分电荷（Gasteiger 电荷）。这两个字段是 PDBQT 的精华，PDB 里根本没有。

4.2.2 原子类型

原子类型决定了 Vina 用什么参数去算范德华力和氢键。常见的几类我列一下：

• C、A、N、O、S —— 对应 sp3/sp2 碳、芳香碳、氮、氧、硫，最常用。
• H、HD、HS —— 氢、极性氢、硫醇氢。别写混了，HD 用来做氢键供体，写错了氢键方向就反了。
• OA、NA、SA —— 带孤对电子的氧、氮、硫，常用于受体端的氢键识别。
• F、Cl、Br、I —— 卤素原子，Vina 1.2 以上版本才支持，低版本会忽略。
• Mg、Ca、Mn、Zn —— 金属离子，需要特殊处理，我一般会额外检查金属配位键。

4.2.3 部分电荷

最后一个数字，通常保留 4 位小数。Vina 用这个电荷计算静电相互作用。如果电荷全为 0，静电项直接废掉，打分结果会少一大块。

标准做法是用 Gasteiger-Marsili 电荷，MGLTools 的 prepare_ligand4.py 会自动计算。我以前试过直接用 MMFF94 电荷往里填，结果对接姿势全变了——老老实实用 Gasteiger 才稳定。

这里给个参考：

• 非极性碳：-0.1 到 +0.1
• 电负性氧：-0.4 到 -0.7
• 氢键供体上的氢：+0.2 到 +0.4

4.3 ROOT / BRANCH / END 分叉树描述符

这部分是 PDBQT 独有的「柔性骨架」。Vina 需要知道配体里哪些键能转、哪些不能转。默认所有键都刚性，但你总不希望对一个氨基酸侧链做对接吧？

4.3.1 ROOT 块

ROOT 标记的是配体的刚性核心。通常是环或不能旋转的骨架段。文件里以 ROOT 开头，ENDROOT 结尾，中间按顺序列出原子。

ROOT
HETATM  1  C1  UNL  1     1.000   2.000   3.000  0.00  0.00   C  0.1000
HETATM  2  C2  UNL  1     1.500   2.500   3.500  0.00  0.00   C  0.0500
ENDROOT

ROOT 块里的原子在对接过程中坐标保持不动，整个配体绕着它们做刚体运动。我个人习惯把苯环、嘧啶这类刚性平面放进 ROOT。

4.3.2 BRANCH 映射

BRANCH 的定义是：BRANCH 起始原子序号 末端原子序号。它告诉 Vina 从哪个键开始可以旋转。每一个 BRANCH 后面跟着对应的原子，直到遇到 ENDBRANCH。

BRANCH   2   3
HETATM  3  O1  UNL  1     1.800   2.800   4.000  0.00  0.00   OA -0.4000
ENDBRANCH

你想想看，Vina 是怎么利用这些信息的？它会在 ROOT 不变的情况下，沿着 BRANCH 的「关节」依次旋转，产生不同的构象。说白了，BRANCH 的数量直接决定了构象搜索的复杂度。

举个例子：

• 如果配体有 0 个可旋转键，Vina 只需要算刚性对接，几十秒搞定。
• 如果配体有 10 个可旋转键，构象空间爆炸式增长，计算时间可能翻几百倍。

4.4 一个完整的 PDBQT 文件长什么样？

咱们看一个真实的配体文件（甲醇类似物），感受一下整体结构：

HETATM 1  C  UNL 1   0.000   0.000   0.000  0.00  0.00   C  0.1000
HETATM 2  O  UNL 1   1.200   0.000   0.000  0.00  0.00   OA -0.4000
HETATM 3  H  UNL 1   1.600   0.800   0.000  0.00  0.00   HD 0.3000
HETATM 4  H  UNL 1  -0.400   0.400   0.800  0.00  0.00   H  0.0000
HETATM 5  H  UNL 1  -0.400   0.400  -0.800  0.00  0.00   H  0.0000
HETATM 6  H  UNL 1  -0.400  -0.800   0.000  0.00  0.00   H  0.0000
ROOT
HETATM 1  C  UNL 1   0.000   0.000   0.000  0.00  0.00   C  0.1000
ENDROOT
BRANCH   1   2
HETATM 2  O  UNL 1   1.200   0.000   0.000  0.00  0.00   OA -0.4000
BRANCH   2   3
HETATM 3  H  UNL 1   1.600   0.800   0.000  0.00  0.00   HD 0.3000
ENDBRANCH
ENDBRANCH
END

注意到没：只有 O-H 键被记为可旋转分支（虽然这里只有 1 个单键旋转，实际上 O-H 键在对接中很少转，Vina 默认把所有非环单键都视为可旋转，你可以通过 torsion 参数手动限制）。

那我再问一个问题：为什么要把 ROOT 单独列出来？因为 Vina 在刚性对接模式下，只优化 ROOT 的位置和朝向；而所有 BRANCH 原子则额外增加内部自由度。分开写，代码逻辑就清晰了。

4.5 知识结构图

下面这张图展示了 PDBQT 文件内部的组织逻辑。建议你打印出来或者贴在常见问题文档里，一看就明白数据从哪里来、到哪里去。

4.6 校验与常见错误

生成完 PDBQT 后，我建议你先跑一遍格式校验。方法很简单：

1. 检查 ROOT 和 ENDROOT 是否成对出现（大小写敏感！）
2. 检查 BRANCH / ENDBRANCH 是否形成平衡的树结构
3. 检查所有原子的类型是否在 AD4_parameters.dat 内
4. 检查部分电荷是否为有效的浮点数（不能写 0.00.00 这种）

我在项目里就用下面的小脚本做快速自检，分享给你参考：

grep -c '^ROOT$' ligand.pdbqt
grep -c '^ENDROOT$' ligand.pdbqt
grep -c '^BRANCH' ligand.pdbqt
grep -c '^ENDBRANCH' ligand.pdbqt

如果 ROOT 和 ENDROOT 数量不相等，或者 BRANCH 和 ENDBRANCH 不匹配，那 Vina 会直接报错退出。别问我是怎么知道的……

4.7 提要与工具推荐

PDBQT 不是你手动去写的格式。我建议用以下工具自动生成：

• MGLTools 的 prepare_ligand4.py：最权威，适合小分子。
• OpenBabel 的 obabel：速度快，支持批量，但类型映射需要调参。
• Meeko：Python 库，可编程，适合批量虚拟筛选流程。

不管用哪个工具，生成后你都应该打开文件看一眼。至少检查一下原子类型是否合理、电荷是否非零、BRANCH 结构是否对称。嗯，在这里多花三分钟，后面对接出错能省三小时调试。