企业官网建设流程全解析

用 AlphaFold 做蛋白质结构比较：RMSD 与 lDDT 的 5 种实用做法

【免费下载链接】alphafoldOpen source code for AlphaFold 2.项目地址: https://gitcode.com/GitHub_Trending/al/alphafold

AlphaFold 是 DeepMind 开源的蛋白质结构预测项目，仓库里除预测流程外，还内置了结构相似度评估的核心代码。围绕蛋白质结构比较这一任务，下面会完成从读取 PDB 坐标、计算 RMSD 与 lDDT，到多模型排名、分区域诊断、批量评估的全部步骤，示例均基于仓库自带模块与测试数据，可直接运行。

什么时候需要比较两个蛋白质结构

跑完一次预测后，或者手上拿到了同一靶点的两份结构（一份来自预测，一份来自实验数据），你通常关心两件事：两者整体差距有多大、哪一份更可信。口头说"看起来差不多"不够有说服力，需要一组可复现的数字。蛋白质结构相似度的量化，基本就靠两个指标：RMSD 和 lDDT。

比较相关的代码集中在两处：alphafold/common/protein.py 负责 PDB 结构的读取与写出，alphafold/model/lddt.py 实现 lDDT 打分。验证代码时可以直接用 alphafold/common/testdata/ 下的2rbg.pdb、5nmu.pdb、glucagon.pdb三个测试文件。

RMSD 与 lDDT：区别与选型

RMSD：叠在灯箱上，看每个地标挪了多远

RMSD（均方根偏差）的思路是"重叠"：先把两个结构叠进同一坐标系，通常要做一次 Kabsch 最优旋转，让整体错位尽可能小；再看对应原子的平均位移。可以把两张透明地图叠在灯箱上理解这件事——无论整体怎么平移旋转，只要每个地标的平均偏移很小，就认为两张地图高度吻合。RMSD 单位是 Å，数值越小，说明两者整体形状越接近。

lDDT：对"距离账本"，而不是对位置

lDDT（local Distance Difference Test）换了参照系：它不关心某个原子在哪，只关心原子之间的相对距离对不对。想象一本"距离账本"，记录了每对残基之间的距离，然后把两本账逐条核对——在 15 Å 的窗口内，某一对距离差值小于 0.5 Å 记满分，0.5~1 Å 记 0.75 分，1~2 Å 记 0.5 分，2~4 Å 记 0.25 分，最终得到 0~1 的加权平均。仓库中 alphafold/model/lddt.py 的lddt()函数按这个分档实现（略去了物理可行性修正项，属于近似版）。由于只看相对距离，整体平移和旋转不会改变分数，因此 lDDT 计算不需要预先做结构对齐。

指标怎么选

一句话决策：想快速判断"是不是同一种折叠"，用 RMSD；想回答"这次预测有多接近真实结构"，用 lDDT。正式评估时建议两个数都记。

动手：安装 AlphaFold 并跑通第一次 lDDT 计算 ⚙️

# 拉取源码 git clone https://gitcode.com/GitHub_Trending/al/alphafold cd alphafold # 安装依赖（含 JAX、Biopython、NumPy 等） pip install -r requirements.txt

下面的最小示例完成三件事：读入两个 PDB、截取公共残基长度、算出第一个 lDDT 分数。

import jax.numpy as jnp from alphafold.common import protein from alphafold.model.lddt import lddt # 读取 PDB 文件，得到 Protein 对象 ref = protein.from_pdb_string(open("alphafold/common/testdata/5nmu.pdb").read()) pred = protein.from_pdb_string(open("alphafold/common/testdata/2rbg.pdb").read()) # atom_positions 形状为 [残基数, 原子类型数, 3] t = jnp.asarray(ref.atom_positions) p = jnp.asarray(pred.atom_positions) n = min(len(t), len(p)) # 示例中截取公共长度 t, p = t[:n][None, ...], p[:n][None, ...] mask = jnp.ones((1, n, 1)) # 标记参考结构中存在的原子 print(float(lddt(p, t, mask))) # 0~1，越接近 1 越接近

注意：lddt要求两个结构按残基一一对应。若序列长度不同或残基顺序不一致，先做序列比对、只取对应残基再计算，否则得到的是没有意义的数字。

任务：给单次预测打分

拿到一份预测结构和一份实验结构后，除了总分，更实用的是定位"哪里不行"。lddt()提供per_residue=True参数，返回每个残基的分数：

per_res = lddt(p, t, mask, per_residue=True)[0] # 每个残基一个分数 overall = float(lddt(p, t, mask)) worst = np.argsort(per_res)[:10] # 得分最低的 10 个残基 print(f"总分 {overall:.3f}，薄弱残基位置 {worst.tolist()}")

解读方式：总分 0.8 以上通常意味着预测质量不错；如果某些残基分数远低于整体均值，多半是 loop 区或侧链取向出了问题，这正是后续放宽或重算时应该关注的区域。

任务：对 5 个候选模型排序

AlphaFold 默认每个靶点输出 5 个模型（入口脚本见 run_alphafold.py）。有参考结构时，按 lDDT 从高到低给模型排序是最常见的筛选方式：

def rank_models(models, ref_pts, mask): """models 为若干预测坐标（与 ref_pts 对齐），按 lDDT 降序返回""" scores = [float(lddt(m[None, ...], ref_pts, mask)) for m in models] return sorted(zip(range(len(models)), scores), key=lambda x: x[1], reverse=True) for idx, s in rank_models(models, t, mask): print(f"model {idx}: {s:.4f}")

没有实验结构可作参考时，lDDT 用不上，此时应参考预测输出自带的 pLDDT 置信度做自评估；lDDT 只负责"两份结构互相对比"这一件事。

任务：批量比较目录下的 PDB 文件

手头有两个目录predictions/与references/，文件同名对应，目标是一次算完所有蛋白质的分数。用多进程并行即可，几十个文件几分钟内跑完：

import os from concurrent.futures import ProcessPoolExecutor def compare_one(name): """比较一对同名 PDB，返回 (文件名, lDDT)""" pred = protein.from_pdb_string(open(f"predictions/{name}").read()) ref = protein.from_pdb_string(open(f"references/{name}").read()) n = min(len(ref.atom_positions), len(pred.atom_positions)) p = jnp.asarray(pred.atom_positions)[:n][None] t = jnp.asarray(ref.atom_positions)[:n][None] return name, float(lddt(p, t, jnp.ones((1, n, 1)))) names = [f for f in os.listdir("predictions") if f.endswith(".pdb")] with ProcessPoolExecutor(max_workers=4) as ex: results = list(ex.map(compare_one, names)) # 4 进程并行

把results存成 CSV 后，后续画图、筛选、汇报都基于这份数据，保证参数一致、结果可复现。

性能与进阶 🚀

分窗口看质量，别只看总分

整条蛋白质的 lDDT 均值会掩盖局部问题。拿到per_residue=True的结果后，按固定窗口切段统计：

W = 40 # 每 40 个残基一个窗口 for s in range(0, len(per_res), W): seg = per_res[s:s + W] print(f"残基 {s}~{s + len(seg) - 1}：{float(seg.mean()):.3f}")

loop 区偏低是常态；如果结构核心区（螺旋、折叠）也明显偏低，说明模型整体有问题，值得换参数重跑。

把批量结果画成分布图

可视化不需要额外工具，一个直方图就能看出整体质量分布和异常点：

import matplotlib.pyplot as plt scores = [s for _, s in results] plt.hist(scores, bins=20, alpha=0.7) plt.axvline(0.8, color="r", ls="--") # 常用参考线 plt.xlabel("lDDT") plt.ylabel("count") plt.show()

参考线位置按项目标准自定，画出来后主要看两件事：分布主体落在哪一段，以及有没有明显偏低的离群点。

计算慢时的三个降速开关

• 只取 Cα 原子：coords[:, 1, :]即每个残基的 Cα 坐标（原子类型索引 1），把 37 种原子类型降到 1 种，距离对数量大幅下降，lDDT 依然有意义。
• 整批传入：lddt支持(batch, length, 3)输入，一次传 5 个模型比循环调 5 次更快。
• 换设备：JAX 会自动选择后端，GPU/TPU 环境下大 batch 的收益明显。

踩坑与排错

• 算出的 lDDT 和评估网站对不上：仓库实现是近似版，lddt.py的注释明确说明省略了物理可行性修正项（键长等约束）；加上 Cα 与全原子的选择差异，不同实现的数值不必要求一致，同一实现下可比即可。
• 两个结构长度不一致：lDDT 的前提是残基一一对应。先做序列比对，只取匹配上的残基，gap 位置直接跳过，不要强行截断。
• RMSD 忽大忽小：几乎可以确定是没做 Kabsch 最优旋转。直接拿原始坐标相减得到的 RMSD 没有任何意义，必须先对齐再算。
• 分数偏低但目测结构很接近：检查参考结构是否本身不完整。lDDT 只统计在参考结构中真实存在的距离对（由 mask 控制），参考数据缺残基、缺原子都会影响分母的构成，压低分数。
• 多链结构无法比较：from_pdb_string支持chain_id参数，先按链拆分，再逐链比较，不要混在一起算。

小结

• RMSD 看全局形状，lDDT 看局部距离，两者口径不同，正式评估建议同时记录。
• 比较成立的前提是"对应关系"：序列先比对、结构先对齐，数字才有意义。
• 从单次打分到批量排名，核心就是from_pdb_string、lddt、per_residue三个入口，其余都是循环、排序和并行。
• 先用 alphafold/common/testdata/ 的测试数据验证整个链路跑通，再接入自己的结构文件。

【免费下载链接】alphafoldOpen source code for AlphaFold 2.项目地址: https://gitcode.com/GitHub_Trending/al/alphafold