📚 课程目录第3章输入文件格式检查

3. 输入文件格式检查：PDBQT规范归纳总结

PDBQT 文件，说白了就是 AutoDock 系列软件的“入场券”。格式对了，任务跑起来才顺；格式错了，要么直接报错，要么跑出来的结果驴唇不对马嘴。我前两年带团队做虚拟筛选的时候，就因为一个配体的ROOT标签写错了位置，整整耽误了两天计算资源——全是白白烧掉的 CPU 时间。后来我养成的第一个习惯，就是所有 PDBQT 文件进队列之前，先拿脚本扫一遍结构。

1. PDBQT 文件的核心结构

一个标准的 PDBQT 文件，由这几大部分组成。我习惯把这个结构叫做“三层夹心”：顶部是注释和参数，中间是原子坐标，底部是分支定义。下面这个表格把每一段的责任写清楚了：

段落	关键词	作用
注释区	`REMARK`	记录来源、原子类型转换信息、电荷参数等
头信息	`ROOT`	声明可旋转体系的刚性部分（配体通常只有一个ROOT）
原子记录	`ATOM` / `HETATM`	每一行的原子坐标、类型、电荷、AD4 type
分支定义	`BRANCH` / `END`	定义可旋转键连接的子树结构
文件结束	`END`	必须独占一行，无额外空格

我记得有一次批处理跑了8个小时，结果发现所有结果的打分都异常高——原来配体文件的ROOT和BRANCH嵌套搞反了，整个柔性被当成了刚体去对接。

2. 原子类型与电荷的规范

PDBQT 里最坑的地方，就是原子类型必须严格匹配 AD4 力场。自动脚本转换的时候，经常把 C 写成 CA 或 CB，区别可大了去了。我给大家列一下核心规矩：

原子类型必须是AD4认识的：C、A、N、OA、NA、HD、Cl、Br、I、P、S、F、Mg、Zn、Ca、Fe…… 自定义类型会直接被AutoDock忽略。
每个原子必须携带 Gasteiger 电荷，写入第55–66列，保留3位小数。
电荷值不能为0.000——至少应该是0.000也算合法，但如果是根本就没有电荷列，比如0.00少了位数，直接崩。

我在项目中遇到过最离谱的一次：蛋白质的Zn原子电荷被转成了2.000，结果对接出来的能量完全不可信。后来我改脚本强制加了个电荷范围检查——超出[-2, 2]直接报错。

3. 可旋转键与柔性残基的处理

ROOT 是整条链的起点。配体通常只有一个 ROOT，但柔性残基可以有多个。写错了的话，AutoDock 很可能读一个原子就报“unexpected token”。下面是一个正确的配体头示例：

REMARK  Name = ligand1
REMARK  AD4 type conversion done by prepare_ligand4.py
ROOT
ATOM      1  C1  LIG     1      11.234   5.678  12.345  0.1234 0.000     C
ATOM      2  O1  LIG     1      12.345   6.789  13.456 -0.5678 0.000     OA
...
ENDROOT
BRANCH   1   3
ATOM      3  C2  LIG     1      13.456   7.890  14.567  0.2345 0.000     C
ENDBRANCH
END

注意：ROOT 和 ENDROOT 必须成对出现，且根原子内部不能包含任何 BRANCH。我曾经见过有人把 ENDROOT 写成了 END ROOT（多了空格），结果整个批次无一例外地解析失败。

4. 常见格式错误与排查方法

搞了这么多年 Autodock，格式错误我闭着眼都能列出来。下面这张表是我自己踩坑的汇总：

错误类型	现象	排查方法
缺少 `ROOT` 或 `ENDROOT`	AutoDock 直接崩溃，报“File is not a valid PDBQT”	用 `grep -c ROOT` 检查标签数量
原子类型不识别	计算出的能量全是 `0.0`	检查 `AD4_type` 列，对照官方类型表
电荷列格式错误	提示“Cannot parse atom record”	用 `awk '{print $6}'` 检查小数点位数
分支嵌套深度超限	AutoDock 只处理前32个可旋转键	用 `grep -c BRANCH` 统计，超过32就得拆分子结构
空白行或异常换行	读取到一半停止，不报错	用 `cat -A` 检查 Windows 换行符（\r\n）

避坑指南：我曾经因为一批 Windows 用户发来的 PDBQT 文件混了 \r\n，用 Linux 的 dos2unix 批量转换后一切正常。

5. 自动化检查脚本思路

手动看几百个 PDBQT 文件不现实，我写了个简单的 bash 脚本，每次提交批次前跑一遍。核心逻辑就三条：数标签、验格式、查电荷。下面是我经常用的版本：

#!/bin/bash
# pdbqt_validator.sh - 快速检查PDBQT格式
for f in *.pdbqt; do
    echo "检查 $f ..."
    # 1. ROOT/ENDROOT必须成对出现
    root_cnt=$(grep -c "^ROOT" "$f")
    endroot_cnt=$(grep -c "^ENDROOT" "$f")
    if [ "$root_cnt" -ne "$endroot_cnt" ]; then
        echo "  [错误] ROOT($root_cnt) 与 ENDROOT($endroot_cnt) 数量不匹配"
    fi
    # 2. 必须有END
    if ! grep -q "^END" "$f"; then
        echo "  [错误] 缺少END行"
    fi
    # 3. 检查电荷列是否为空或非数字
    awk '/^ATOM|^HETATM/ {chg=substr($0,55,12); if(chg+0==0 && chg!="0.000") print FILENAME": 电荷异常 "chg}' "$f"
done

这个脚本并不能查出所有问题，但能把绝大多数导致批次崩溃的“硬伤”揪出来。我每次跑完筛选项目都会顺手更新一下规则，现在已经有20多个检查项了。

嗯，格式检查这事儿，说难不难，但一旦疏忽，代价就是好几个计算天的浪费。把上面的结构记牢，再配个自动化脚本，基本能杜绝因为 PDBQT 格式造成的批次质控失败。