📚 课程目录第3章 PDB文件解析与处理

PDB文件解析与处理

PDB文件是结构生物学的通用语言。你拿到的结构数据，不管是X射线、核磁还是冷冻电镜，最终大多以PDB格式存下来。但PDB文件说白了就是个文本表格，一行行记录原子坐标、残基信息、温度因子……所以我们需要一个趁手的工具来解析它。我个人几乎每天和Bio.PDB打交道，今天就聊聊怎么用Biopython读取、写入和操作PDB文件。

为什么是Biopython？

Biopython的Bio.PDB模块是专门干这活的。它把PDB文件解析成一个层级树：Structure → Model → Chain → Residue → Atom。每个节点都有属性和方法。我记得刚开始学时总忘这些层级，后来写了个小脚本打印层级结构，才彻底搞清楚。

读取PDB文件

读取是最基础的操作。用PDBParser，一行代码就能把一个PDB文件转成Structure对象。我通常这么写：

from Bio.PDB import PDBParser
parser = PDBParser(QUIET=True)  # QUIET抑制警告
structure = parser.get_structure("1ubq", "1ubq.pdb")

这里QUIET=True很重要。如果PDB文件有一些格式瑕疵（比如非标准残基名），解析器会抛警告。设成True省得屏幕刷屏。我在处理批量PDB时吃过亏——没设QUIET，结果日志文件好几GB。

拿到structure后，你可以遍历所有model、chain、residue、atom。例如打印所有链ID：

for model in structure:
    for chain in model:
        print(chain.get_id())

小技巧：如果PDB文件包含多个模型（比如NMR结构），structure会包含多个Model对象。用structure[0]获取第一个模型。

写入PDB文件

修改完结构，得写回PDB文件。Bio.PDB提供了PDBIO类。用法也简单：

from Bio.PDB import PDBIO
io = PDBIO()
io.set_structure(structure)
io.save("output.pdb")

不过默认写入会保留所有原始ATOM记录。如果你只想写某个链或残基，可以先提取子结构再保存。我记得有一个项目需要把多链复合物拆成单个文件，我就是用io.set_structure(chain)然后保存。

操作PDB结构

操作包括添加、删除、修改原子和残基。以删除一个残基为例：

# 假设要删除第一个模型的A链的第10号残基
model = structure[0]
chain = model['A']
residue = chain[(' ', 10, ' ')]  # (hetero, seqnum, insertion)
chain.detach_child(residue.get_id())

注意残基的标识符是一个元组：(异标, 序列号, 插入码)。比如水分子是('W', 1, ' ')，标准残基是(' ', 10, ' ')。这里我踩过坑：一开始我直接用整数索引，结果抛KeyError。所以查残基一定得用这个三元组。

修改原子坐标更常见。比如平移整个结构：

from Bio.PDB import Vector
translation = Vector(5.0, 0.0, 0.0)
for atom in structure.get_atoms():
    atom.set_coord(atom.get_vector() + translation)

注意：坐标操作后记得重新保存。但不要直接用同一个文件覆盖原PDB，以免丢失原始数据。我习惯在文件名加个_backup后缀。

常用属性提取

做结构分析时经常需要提取B因子、残基名称、原子名称等。下面是我常用的几个代码片段：

# 获取所有原子坐标和B因子
coords = []
bfactors = []
for atom in structure.get_atoms():
    coords.append(atom.get_vector())
    bfactors.append(atom.get_bfactor())

# 获取残基序列（单字母代码）
from Bio.PDB.Polypeptide import PPBuilder
ppb = PPBuilder()
for pp in ppb.build_peptides(structure):
    print(pp.get_sequence())

PPBuilder可以自动从结构中提取多肽链的序列。我经常用它和序列比对工具配合，检查PDB中的序列与实验序列是否一致。

处理缺失原子与交替构象

晶体结构中常有缺失的loop或侧链。Biopython的Atom对象可能缺失。你可以用atom.is_disordered()判断交替构象。我曾经用以下方式获取占据率最高的构象：

if atom.is_disordered():
    selected = atom.selected_child
else:
    selected = atom

  重点：PDB解析时默认会保留所有交替构象。如果你只想要主链，可以先做选择。用Polypeptide.PPBuilder只构建主链多肽，不需要关心侧链。

常见操作速查表

操作	代码示例	说明
解析PDB	`parser = PDBParser(QUIET=True) structure = parser.get_structure("id", "file.pdb")`	QUIET=True 抑制格式警告
遍历链	`for chain in model: print(chain.get_id())`	model是Structure的第一个模型
获取残基	`res = chain[(' ', 10, ' ')]`	标识符为(异标, 序列号, 插入码)
删除残基	`chain.detach_child(res.get_id())`	从链中移除该残基
修改坐标	`atom.set_coord(new_vector)`	new_vector是Vector对象
写入PDB	`io = PDBIO(); io.set_structure(s); io.save("out.pdb")`	可保存子结构（链或模型）

嗯，基本就这些。Biopython的PDB模块功能很全，但常用操作就这么几板斧。你上手之后，会发现它比手动解析PDB文本快得多。我始终觉得，学会解析和操作PDB是做结构生物信息学分析的第一步——后面做对齐、做相似性比较都依赖这个基础。