📚 课程目录第3章输出文件解析

输出文件解析：对接结果文件

跑完Vina，你会拿到两个关键文件：一个pdbqt，一个log。很多人只看log里的打分，觉得构象文件无所谓——我以前也这么干过，后来发现这样会漏掉不少信息。今天就掰开揉碎，把这两个文件里的门道说清楚。

PDBQT 结果文件：构象的“全家福”

Vina每次对接都会输出一个多模型PDBQT文件。说白了，就是把所有找到的构象塞进一个文件里，一个接一个排好。我习惯叫它“全家福”，因为这张照片里有多少个成员，完全取决于你设的--num_modes参数。

字段逐行解析

打开一个典型的pdbqt，你会看到类似这样的结构：

MODEL 1
REMARK  Name = ligand_01.pdbqt
REMARK  VINA RESULT: -9.5 0.0 0.0
ATOM      1  C   LIG     1      12.345   56.789   23.456  0.00  0.00    C
ATOM      2  C   LIG     1      12.567   56.901   23.234  0.00  0.00    C
...
ENDMDL
MODEL 2
REMARK  VINA RESULT: -9.2 1.2 2.3
...
ENDMDL

我来拆解一下关键行：

MODEL 1 / ENDMDL：构象编号。每个模型代表一个独立的对接姿态。编号从1开始，最大就是你设的num_modes。不过有一点要注意——Vina可能会少输出几个，如果它觉得剩下的构象能量太高或者太相似。
REMARK VINA RESULT：这行是灵魂。三个数字依次是：对接打分（kcal/mol）、与最优构象的RMSD、与最优构象的上半部分RMSD。我在项目中经常用第二个数值来判断“这个构象和最好的那个有多像”。
ATOM 行：和输入pdbqt的原子描述格式完全一样，只是坐标变了——现在这个坐标是对接后的位置。注意最后两列：第一列是电荷，第二列是AD4类型。如果是0.00，说明这个原子在对接过程中没有被考虑vdW或静电作用。

避坑指南：我曾经遇到过一个案例——对接大分子时，MODEL行前面多了空格，导致某些脚本按行解析时漏掉了第一个模型。建议写解析脚本时用正则匹配 ^MODEL\s+\d+，别偷懒用简单的split。

怎么提取单个构象？

两个主流办法：

手动分离：用文本编辑器搜索 MODEL 3，把从这行到对应的 ENDMDL 之间的内容复制出来，另存为一个单独的pdbqt。适合临时看一两个构象。
脚本批量拆分：我常用一个简单的Python脚本：

import re
with open('output.pdbqt', 'r') as f:
    content = f.read()
models = re.split(r'ENDMDL\n?', content)
for i, model in enumerate(models, 1):
    if model.strip():
        with open(f'model_{i}.pdbqt', 'w') as out:
            out.write(model.strip() + '\nENDMDL\n')

脚本会自动把每个构象存成单独文件。命名带上编号，方便后续批量处理。

Log 文件：对接的“体检报告”

Log文件记录了对接过程的详细信息。它不是必需的输出来，但很多细节只有这里能读到。这是它的标准格式：

mode |   affinity | dist from best mode
     | (kcal/mol) | rmsd l.b.| rmsd u.b.
-----+------------+----------+----------
   1        -9.5      0.000      0.000
   2        -9.2      1.234      2.345
   3        -8.7      2.456      3.678
   4        -8.5      2.890      4.012
   5        -8.3      3.456      4.789
   6        -8.1      3.789      5.012
   7        -7.9      4.234      5.678
   8        -7.5      4.890      6.012
   9        -7.3      5.456      6.789
  10        -7.0      5.890      7.012

表格的列很好理解：

mode：构象编号，和pdbqt里的MODEL编号严格对应。第1号就是打分最好的那个。
affinity：预测的结合自由能。越负越好。我一般觉得-8以下可以进入候选名单，-10以上算比较优秀。
rmsd l.b.：与最佳构象的RMSD（Lower Bound）。这个值越小，说明两个构象越接近。如果两个构象的RMSD小于1.0Å，我通常认为它们算同一类。
rmsd u.b.：与最佳构象的RMSD（Upper Bound）。通常比l.b.大一点，提供一个更宽松的比较范围。

个人经验：我习惯把log文件的第一列（mode）和第三列（rmsd l.b.）贴到Excel里，快速筛选出“打分靠前但构象差异明显”的候选。这样既能保证能量好，又能拿到多样性的构象——筛选命中率的常用策略。

Log里面的隐藏信息

除了那个表格，log文件前面还有几段，很多人会跳过。我建议你养成看的习惯：

种子信息：记录了你使用的随机种子。如果后续需要复现结果，这个值必须记录下来。
搜索空间尺寸：会告诉你box的中心和大小。我曾经有一次改过结构后忘了更新box，结果log里的坐标一看还是旧的……幸好日志帮我定位了问题。
对接跑了几轮：Vina默认是8轮穷举搜索。如果设置exhaustiveness=4，这里会显示4轮。

构象提取与排序：从数据到“干活”

拿到构象和打分后，怎么选？这是我的习惯：

看打分：优先选top 3的构象。它们能量最低，理论上最可能发生。
看构象差异：如果top 3的RMSD都小于1Å，说明它们几乎一样，那就只提取一个。如果差异很大（RMSD>2Å），就把所有明显不同的构象都提取出来。
做聚类分析：这是高级玩法。我会对10个构象做RMSD聚类，分成3-4个集群，每个集群选一个代表。这样做出来的结果既全面又互斥。

之所以强调这一步，是因为你对接的目的是找到真实的结合模式，不是找一个“分数最高的构象”。我自己就踩过坑——有一次只提取了top 1构象，结果后面做MD模拟发现它很快就跑了，反而是排名第3的构象在动力学中保持稳定。

核心结论：pdbqt文件存的是构象坐标，log文件存的是评分和相似度。两个要结合看，缺一不可。提取构象时，别只盯着分数，多样性往往是发现真结合模式的关键。

你可能会问：那log里的RMSD到底用哪个？我一般只看RMSD l.b.，它对公共骨架变化敏感度适中，能快速筛选出明显不同的构象。如果l.b.超过2Å，说明结合模式大概率不一样。至于10个构象还不够用的情况——嗯，那就加大num_modes，不过要记住，Vina默认只输出9个，最多能设到20。

总之，读懂这两个文件是Vina进阶的第一步。下次遇到对接结果看起来很诡异时，别急着改参数，先回去看看log里有没有异常值，或者pdbqt里是不是混入了奇怪的东西。我宁愿多花5分钟检查，也好过跑完一轮才发现问题出在数据理解上。