企业官网建设流程全解析

如何判断 AlphaFold 预测的蛋白质结构靠不靠谱？RMSD 与 lDDT 结构相似度实用指南

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

拿到 AlphaFold 的预测 PDB 文件后，下一个问题总是：它跟实验结构差多远？肉眼叠合能看个大概，但汇报结论需要数字。这篇文章从结构比较的对齐步骤讲起，用可运行的代码实现 RMSD 与 lDDT 两个最常用的相似度指标，最后给一份选型与避坑清单。

比较之前先对齐：原子配对、旋转对齐与原子子集

直接对两个 PDB 文件做坐标相减没有意义——同一个蛋白质在坐标系里可以放在任何位置、转任意角度。真正能做比较，要先完成三件事：

原子配对。先决定预测结构里的哪个原子对应实验结构里的哪个原子。序列完全一致时按残基序号一一对应；有插入缺失时，先跑一次序列比对，只保留共同区域。

平移与旋转对齐。把两个结构的质心都挪到原点，再找让总距离最小的旋转。这一步一般用 Kabsch 算法（奇异值分解 SVD）解出最优旋转矩阵，也就是各种结构软件里 superpose 按钮背后的东西。

原子子集选择。全原子比较会被侧链的细碎抖动主导，工程上常只对 Cα 原子或主链打分。AlphaFold 内部用一张固定 37 格表存放每个残基的原子（见 residue_constants.py 中的atom_order，CA 在 1 号位），"取 Cα 那一列"只是一个索引的事。

图里分叉的两条路，就是后两章的两把标尺：RMSD 走"先对齐再相减"，lDDT 走"不对齐只比距离"。

全局标尺：从 RMSD 公式到 Kabsch 对齐代码

一句话：均方根偏差 RMSD（Root Mean Square Deviation）就是两个结构做最优对齐后，每对对应原子的"平均偏离"。

打个比方：它像比对两张印在同一城市上的纸质地图——先把其中一张平移、旋转，让地标尽量重合，再算每个地标偏离多少。数值越大，说明两张地图越"对不上"。

数学形式是平方距离平均后开根号：

RMSD = sqrt( (1/N) · Σᵢ ‖xᵢ − x'ᵢ‖² )

工程实现里唯一有技术含量的地方，是求"最优旋转"。下面是仅依赖 NumPy 的完整实现：

import numpy as np def kabsch_rmsd(pred: np.ndarray, target: np.ndarray) -> float: """Kabsch 最优旋转对齐后计算 RMSD；pred/target 形状为 [N, 3]""" # 1) 质心归零：两个结构都平移到原点 pred_c = pred - pred.mean(axis=0) targ_c = target - target.mean(axis=0) # 2) 交叉协方差矩阵做 SVD（行向量约定） u, _, vh = np.linalg.svd(pred_c @ targ_c.T) v = vh.T # 3) 防"镜像"：det<0 说明解出了反演旋转，翻转最后一主轴 d = np.diag([np.linalg.det(v @ u.T), 1.0, 1.0]) rot = u @ d @ v.T # 4) 应用旋转，算平方距离的均方根 delta = pred_c @ rot - targ_c return float(np.sqrt(np.mean(np.sum(delta ** 2, axis=1))))

rot = u @ d @ v.T是 SVD 分解出来的闭式解；d处理行列式为负的镜像情形，避免把某个结构"翻"成手性相反的镜像。

为什么工程上常只对 Cα 计算？主链 Cα 链就像项链的"串绳"：它决定整体折叠，而原子数量只有全原子的零头。侧链数量多、摆动大，混进来会把指标搅浑。顺带一提，仓库里结构弛豫步骤（relax.py）用 RMSD 的方式更简单：直接相减最小化前后的两组坐标（源码第 70 行），观察结构在能量精修中"动了多少"。因为比较的是同一套坐标的两个时刻，连对齐都省了。

局部标尺：lDDT 免对齐打分原理与实现

RMSD 有三类失效场景：

• 局部变化被平均掉：某个环区打开 1 Å，摊到整条链上几乎看不见，"平均"天然掩盖局部问题；
• 原子缺失：实验结构里有些残基没解析出坐标，预测结构却是满的，直接相减数量都对不上；
• 序列长度不同：同源蛋白长度不一致，无法原子对原子地比。

局部距离差异测试 lDDT（local Distance Difference Test）把这三个问题全绕开了：它完全不看绝对坐标，只看距离。

类比一下：不必把两张地图对齐，只要逐条核对"体育场离市政厅 5 公里"这样的局部距离关系是否一致。所有局部距离都对得上，局部结构就是可靠的。所以 lDDT 也常被称为"superposition-free（免叠合）"指标。

核心思路四步：① 分别构建真实与预测结构的距离矩阵；② 圈出真实结构中距离小于 15 Å 的残基对，屏蔽缺失区域与自身对；③ 对每一对算 L1 距离差；④ 按固定分箱（<0.5 / 1 / 2 / 4 Å 各计 0.25 分）打分后归一化。仓库的 lDDT 源码 就是这套逻辑的 JAX 实现，lddt_test.py 的单元测试逐一验证了分箱行为：偏差 0.5 Å 内满分，超过 4 Å 就归零。

下面是与源码一一对应的 NumPy 版：

import numpy as np def lddt_score(pred: np.ndarray, target: np.ndarray, mask: np.ndarray, cutoff: float = 15.0) -> np.ndarray: """简化版 lDDT；pred/target [batch, N, 3]，mask [batch, N, 1]（1 表示原子存在）""" # 1) 真实 / 预测两个结构的距离矩阵（加 1e-10 防止除零） dt = np.sqrt(np.sum((target[:, :, None] - target[:, None, :]) ** 2, -1) + 1e-10) dp = np.sqrt(np.sum((pred[:, :, None] - pred[:, None, :]) ** 2, -1) + 1e-10) # 2) 只给真实结构中距离小于 cutoff、双方都存在且非自身的残基对打分 eye = np.eye(dt.shape[1], dtype=np.float32) score_pairs = ((dt < cutoff) * mask * mask.transpose(0, 2, 1) * (1.0 - eye)) # 3) 两个结构同对原子的 L1 距离差 l1 = np.abs(dt - dp) # 4) 四个固定分箱 <0.5 / 1 / 2 / 4 Å，每箱各贡献 0.25 s = 0.25 * ((l1 < 0.5) + (l1 < 1.0) + (l1 < 2.0) + (l1 < 4.0)) # 5) 对参与打分的点对归一化 return (score_pairs * s).sum(axis=(-2, -1)) / (score_pairs.sum(axis=(-2, -1)) + 1e-10)

score_pairs一行是全函数的心脏：它决定"哪些点参评"（截断 + 掩码 + 排除自作用）。0.25 * (...)则是分箱打分本身。

🧭 RMSD 与 lDDT 选型指南：多模型批量评估实战

场景化决策，每句一行：

• 比较同一蛋白的两个模型、判断整体是否正确 → Cα RMSD，快且直观；
• 没有实验结构、想评估预测本身质量 → 看 pLDDT，即"预测出来的 lDDT"，定义见 confidence.py 与仓库技术报告；
• 结构有缺失区域、序列长度不一致 → lDDT，RMSD 根本算不出来；
• 怀疑某个结构域塌了 → 逐残基 lDDT，能定位到具体残基；
• 批量筛选上千个模型 → 先 Cα RMSD 粗筛，存活者再算 lDDT 精排。

仓库自带的预测 vs 实验对照如下，每个目标都标了 GDT 分数（与 lDDT 同属距离矩阵家族）：

实际工作里经常要把多个预测模型批量跟一个实验结构比，并同时给出两个指标：

import numpy as np import matplotlib.pyplot as plt def evaluate(models, target, ca_index=1): """批量比较多个预测模型与实验结构，返回 (rmsd_list, lddt_list)""" rmsds, lddts = [], [] for m in models: p, t = m[:, ca_index, :], target[:, ca_index, :] # 只取 Cα rmsds.append(kabsch_rmsd(p, t)) # 全局指标 mask = np.ones((1, t.shape[0], 1), dtype=np.float32) lddts.append(float(lddt_score(p[None], t[None], mask))) # 局部指标 return rmsds, lddts # 模拟 5 个模型：真实结构叠加不同强度的噪声 rng = np.random.default_rng(42) target = rng.normal(size=(80, 37, 3)).cumsum(axis=0) / 8 + 10 models = [target + rng.normal(scale=s, size=target.shape) for s in (0.1, 0.4, 0.8, 1.5, 3.0)] rmsds, lddts = evaluate(models, target) fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4)) ax1.bar(range(1, 6), rmsds); ax1.set_title('CA RMSD / Å（越小越好）') ax2.bar(range(1, 6), lddts); ax2.set_title('lDDT（越大越好）') plt.tight_layout(); plt.savefig('model_comparison.png', dpi=120)

evaluate里对每个模型各调一次 Kabsch 与 lDDT，两个指标共用同一份 Cα 抽取。画出来的两条曲线讲的是同一个故事，但"语气"不同：RMSD 从 0.1 Å 噪声就开始上涨，lDDT 在 0.5 Å 内仍近乎满分。小偏差看 RMSD 更灵敏，大偏差与局部崩坏看 lDDT 更直接；两者打架时，用逐残基 lDDT 定位差异区域。

⚠️ 踩坑清单：结构比较的常见工程问题

📌 要点速查

• 算 RMSD 先做 Kabsch 对齐；对齐做不了或序列长度不同，选 lDDT。
• 快速看全局：Cα RMSD；看局部质量、带缺失原子：lDDT（逐残基）。
• 先建掩码再计算，缺失原子是脏数据第一大来源。
• 汇报 RMSD 时写清楚用了哪些原子、是否已对齐，否则数字不可比。
• 没有实验结构时，pLDDT 是 lDDT 思想的预测版，看分布区间而不是单点均值。

结构比较的下一站，是把标尺从单帧静态快照搬到轨迹层面——动态构象下的相似度指标，大概会延续距离矩阵这条路线继续走远。

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