4. 蛋白结构数据获取与预处理归纳总结

蛋白结构数据，是药物发现的起点。

数据质量不行，后面模型再强也没用。我见过太多项目，数据准备不到位，训练出来的模型根本没法用。说白了，这一步决定了整个项目的上限。

4.1 蛋白结构数据的主要来源

目前最主流的蛋白结构数据库有三个。我个人习惯按照使用场景来选择：

数据源/格式	特点说明
PDB (Protein Data Bank)	实验测定结构，权威但覆盖率有限。PDB中约90%是X-ray，其余为NMR、Cryo-EM和少量计算模型。适合有明确实验结构的靶点。
AlphaFold DB	DeepMind预测的高质量结构库，覆盖蛋白质组尺度。最近一次更新后，预测精度在多个关键指标上已接近实验结构。适合缺乏实验数据的靶点。
UniProt + PDB匹配	通过序列ID或同源性匹配获取结构，常用于批量数据准备。需要注意序列一致性阈值，一般建议>90%才能直接使用。

4.2 数据获取的操作流程

以PDB为例，获取结构数据通常有两种方式。我建议初学者先手动试试，熟悉数据格式后再上编程接口。

1. 手动下载：访问 rcsb.org，输入PDB ID（如4AG8），选择PDB格式文件。
2. 程序化获取：使用BioPython或Pymol脚本批量下载。
3. AlphaFold批量下载：利用UniProt ID拼接URL直接获取。

下面给一个批量下载的Python示例。代码不多，但很实用：

import urllib.request
import os

# PDB批量下载
pdb_ids = ['4AG8', '4AG7', '5T91']
base_url = 'https://files.rcsb.org/download/'

for pid in pdb_ids:
    url = base_url + pid + '.pdb'
    urllib.request.urlretrieve(url, f'./pdb_files/{pid}.pdb')
    print(f'Downloaded {pid}')

# AlphaFold批量下载
uniprot_ids = ['P12345', 'Q6GZX4']
af_url_template = 'https://alphafold.ebi.ac.uk/files/AF-{}-F1-model_v4.pdb'

for uid in uniprot_ids:
    url = af_url_template.format(uid)
    urllib.request.urlretrieve(url, f'./af_files/{uid}.pdb')

4.3 预处理的核心步骤

拿到原始PDB文件后，不能直接用。为什么？原始文件中包含大量杂质：水分子、配体、辅因子、溶剂分子，还有可能链标识不完整。

预处理管线（Pipeline）通常包含以下步骤：

• 去除水分子：HOH记录行直接删除。
• 保留配体/辅因子：依据配体3字母代码（如ATP、FAD）保留，但非药物相关小分子（如DMS）建议去除。
• 链的选择：只保留与配体结合相关的链。
• 氢原子添加：需使用工具（如PDBfixer、Reduce）补全，用于静电计算。
• 残基编号重排：部分PDB文件残基编号不连续，需要标准化。

举个具体的操作例子，我曾经就在一个项目里遇到PDB文件链标识缺失，导致后续识别失败。后来我养成了每次预处理都检查链标识的习惯。

使用PDBfixer进行氢原子添加和结构调整：

# 通过命令行调用PDBfixer
pdb_fixer my_protein.pdb --add-hydrogen --output fixed_protein.pdb

# 去除水分子和特定配体
pdb_fixer my_protein.pdb \
  --remove-waters \
  --remove-ligand DMS \
  --output cleaned_protein.pdb

4.4 数据验证与质量控制

预处理完毕，不代表数据合格。你想想看，数据质量差，模型学到的就是噪声。

验证指标和标准：

检查项	验收标准	工具/方法
残基覆盖率	序列覆盖度 > 90%	BLAST比对
立体化学质量	Ramachandran favored > 98%	MolProbity/ProCheck
溶剂可及表面积	无异常暴露的疏水残基	FreeSASA/Areaimol
结构完整度	无非标准残基、断链	PDB格式检查器
分子间冲突	clash score < 5%	MolProbity

我建议写成脚本自动跑一遍。不然几百个结构一个个看，眼睛会瞎掉。

import procheck
from Bio.PDB import PDBParser

parser = PDBParser()
structure = parser.get_structure('protein', 'cleaned_protein.pdb')

# 运行ProCheck
report = procheck.run_procheck(structure)
print(f'Ramachandran favored: {report.favored_percent:.2f}%')
print(f'Clash score: {report.clash_score:.2f}')

# 完整性检查
model = structure[0]
residues_found = len([r for r in model.get_residues()])
expected_residues = 350  # 根据UniProt序列计算

if residues_found / expected_residues < 0.9:
    print('Warning: Low coverage!')
else:
    print('Coverage check passed.')

4.5 特征提取的数据准备

预处理完成的结构，需要转换为特征提取的输入格式。最常见的做法是生成结构描述符文件或直接提取坐标。

常用格式转换：

• 结构描述符：如MACCS指纹、ECFP4（针对蛋白表面性质），由PyBioMed或ProLIF生成。
• 坐标向量：提取Cα原子或侧链重原子的3D坐标，用于图神经网络。
• 距离矩阵：计算残基间Cα距离，生成接触图（contact map）。接触图可用于蛋白质-配体相互作用预测。

接触图的生成代码示例：

import numpy as np
from Bio.PDB import PDBParser

def calc_contact_map(pdb_file, cutoff=8.0):
    parser = PDBParser()
    structure = parser.get_structure('protein', pdb_file)
    model = structure[0]

    # 获取所有Cα原子坐标
    coords = []
    for chain in model:
        for residue in chain:
            if residue.has_id('CA'):
                coords.append(residue['CA'].get_vector().get_array())
    coords = np.array(coords)

    # 计算距离矩阵
    n_res = len(coords)
    dist_matrix = np.linalg.norm(coords[:, np.newaxis] - coords[np.newaxis, :], axis=-1)

    # 生成接触图
    contact_map = (dist_matrix < cutoff).astype(int)
    return contact_map

# 示例使用
contact = calc_contact_map('cleaned_protein.pdb', cutoff=8.0)
np.save('contact_map.npy', contact)
print(f'Contact map shape: {contact.shape}')

这一步的输出，就是喂给机器学习模型的“料”。

嗯，到这里，整个数据获取与预处理的管线就完整了。我把核心逻辑画成了一张流程图，方便你整体把握。

整个流程走下来，其实就是：数据源选择→批量获取→预处理（去水、选链、加氢）→质量验证→特征格式输出。

记住，数据准备阶段投入1小时，后面建模能省10小时。我在这上面吃过亏，现在不管多忙都严格按流程走。