十年匠心定制 · 商业建站与技术教学双线并行 咨询热线:400-886-1026 service@lmnt.cn
ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

蛋白质结构比较:RMSD和lDDT怎么选,才不把相似性算错

蛋白质结构比较:RMSD和lDDT怎么选,才不把相似性算错 蛋白质结构比较RMSD和lDDT怎么选才不把相似性算错【免费下载链接】alphafoldOpen source code for AlphaFold 2.项目地址: https://gitcode.com/GitHub_Trending/al/alphafoldAlphaFold 是一个开源的蛋白质结构预测项目给它一段序列它给你一套三维原子坐标。拿到预测结果后第一个问题几乎总是预测得准不准——但准本身就是个需要量化的概念。用实验结构叠在预测结构上看两者看起来完全一样可一算分数RMSD 给到 1.4 ÅlDDT 只有 0.71。这不是谁的 bug而是两把尺子量的东西不一样一把看整体重合一把看局部接触。本文用 AlphaFold 的 lDDT 源码 和几段可以直接跑的代码把蛋白质结构比较里最核心的 RMSD、lDDT、Kabsch 对齐讲清楚。看起来一样却分数对不上两把尺子的分工先给结论RMSD 量的是叠在一起后每个原子离了多远lDDT 量的是结构内部哪些点挨得近这件事有没有被破坏两者不矛盾测的维度不同。把蛋白质想成一串乐高积木。RMSD 的玩法把两套积木叠放到最重合的位置然后逐块测量同一块积木中心点之间的位移最后求均方根。它关心的是对应原子第 37 号残基对第 37 号残基错位了就如实记账。lDDT 的玩法不叠放。分别画出两套积木的邻居地图——谁和谁挨得近比如 15 Å 以内的原子对然后逐条核对真实结构里这俩挨着预测结构里还挨着吗挨着的程度差了多少它关心的是内部距离关系所以天然不需要对齐两套坐标。同一张折叠图末端螺旋整体转了 20°RMSD 会被这几个对应原子的位移拉高因为它是全局记账而 lDDT 只看到末端局部接触略有变化整体分数掉得少。反过来核心区域两个残基悄悄换了相对位置RMSD 可能只动一点点lDDT 却立刻记上一笔。看到两套分数打架先别怀疑计算先想清楚你关心的是全局形状还是局部接触。RMSD 的 4 步计算流程结论RMSD 本身不含对齐对齐是它的前置步骤而 Kabsch 对齐求让两个点云重合的最优旋转的方法就是最标准的做法。选原子常用 Cα每个氨基酸残基的主链 α 碳既能代表骨架又省计算。平移到质心两套坐标各自减去质心消除位置差异。Kabsch 旋转用奇异值分解SVD把矩阵拆成旋转×缩放×转置旋转的方法求最优旋转矩阵 R让两套点云尽量重合。求均方根对应原子位移差的均方根。最简公式$$RMSD \sqrt{\frac{1}{N}\sum_{i1}^{N} d_i^2}$$这句话的意思是N 对原子每对算一个位移 d_i平方、平均、开根号——平方让大的错位被重罚开根号把量纲拉回 Å。lDDT 评分规则详解距离矩阵与 4 档打分结论lDDT 不碰对齐它把两套结构各自展开成距离矩阵两两原子间距离组成的表格只核对真实结构里距离小于 cutoff 的原子对按偏差分 4 档记分。最简公式$$lDDT \frac{1}{N_{cut}}\sum_{i,j cutoff} \mathbb{1}\big(|d_{ij}^{pred} - d_{ij}^{true}| 4\ \text{Å}\big)$$逐符号解释d^true_ij / d^pred_ij真实结构和预测结构里第 i、j 两个原子各自的距离N_cut真实结构里距离小于 cutoffAlphaFold 取 15 Å的原子对总数只数得上的对才参与评分|d^pred − d^true|预测把这段局部距离改掉了多少15 Å cutoff只关心局部接触两条远端的链再怎么挪都不计分4 档记分真实 lDDT 用 0.5/1/2/4 Å 四档偏差越小分越高所以单对满分是 0.25×41。上面的公式只写了 4 Å 这一档方便阅读完整版见下一节源码。源码精读AlphaFold 的 lddt 函数在为什么设计买单实现只有 alphafold/model/lddt.py 不到 90 行三个设计点值得停下来想。1cutoff 只按真实结构判定再排除自身、乘上掩码。评分范围由真实结构里谁和谁挨着决定而不是预测结构——否则一个整体塌掉的预测会凭空制造一堆新接触反而多得分。掩码 true_points_mask 则负责把缺失原子实验结构里没解析出来的区域从评分里剔除。2为什么处处加 1e-10归一化要除以有效原子对数某个残基若恰好没有落在 cutoff 内的接触分母是 0。加 epsilon 让它在 JAX 里可微、不产生 NaN。3per_residue 为什么总分校准不到每个残基每个残基能凑出的接触对不同核心残基多、末端残基少所以总分数不等于逐残基分数的简单平均。这也是 pLDDT 曲线要按残基画、而不是只看一个总数原因。模型调用侧alphafold/model/modules.py 第 1093 行能看出真实用法取第 1 个原子Cα、带 all_atom_mask、per_residueTrue分数还会被分箱后变成网络头部的训练目标# 训练时 lDDT 的分箱目标modules.py 简化 bin_index jnp.floor(lddt_ca * num_bins).astype(jnp.int32) bin_index jnp.minimum(bin_index, num_bins - 1) # lddt_ca1 时防越界 lddt_ca_one_hot jax.nn.one_hot(bin_index, num_classesnum_bins)动手实践用 NumPy 重算两套分数不装 JAX 也能验证核心逻辑下面两个函数约 30 行覆盖 Kabsch 对齐和 4 档记分可直接复制运行。import numpy as np def ca_rmsd(a, b): a, b a - a.mean(0), b - b.mean(0) # 1) 平移到质心 u, _, vt np.linalg.svd(a.T b) # 2) Kabsch: SVD 求最优旋转 d np.sign(np.linalg.det(u vt)) # 3) 修正手性确保是旋转而非反射 R u np.diag([1, 1, d]) vt return float(np.sqrt(np.mean(np.sum((a R - b) ** 2, 1)))) # 4) 均方根 def mini_lddt(true_pts, pred_pts, cutoff15.): D np.sum((np.asarray(true_pts)[:, None] - np.asarray(true_pts)[None]) ** 2, -1) ** .5 Dp np.sum((np.asarray(pred_pts)[:, None] - np.asarray(pred_pts)[None]) ** 2, -1) ** .5 pairs (D cutoff) ~np.eye(len(D), dtypebool) # 只评 cutoff 内的真实接触 diff np.abs(D - Dp) bins (diff .5) (diff 1.) (diff 2.) (diff 4.) # 4 档记分 return float((0.25 * bins)[pairs].mean())输入输出两个[N, 3]坐标数组输出一个 float。用 alphafold/model/lddt_test.py 里的参数做交叉验证2 个原子、真实距离 5 Å预测距离分别挪到 5.5/6/7/9/11 Å结果依次是1.0 / 0.75 / 0.5 / 0.25 / ≈0——正好每差一档掉 0.25评分逻辑立刻可检验。典型报错与处理AssertionError: len(predicted_points.shape) 3直接调lddt.lddt()时输入必须是(batch, length, 3)记得pts[None, ...]补一个 batch 维ValueError: 形状不匹配Cα 是每残基 1 个点不是 37 个原子——取错索引Cα 是第 1 个原子见 alphafold/common/residue_constants.py 的 atom_order就会撞上长度断言。RMSD vs lDDT 选型对照我该用哪个结论先问自己三个问题答案基本就定了——结构完整吗比的是全局形状还是局部细节数据里有没有缺失原子维度RMSDlDDT量的是什么对齐后对应原子的位移结构内部距离关系需要 Kabsch 对齐是否superposition-free取值≥0 Å越小越好0~1越大越好对末端/局部扰动敏感全局记账局部记账整体更稳缺失原子要手动选同一组原子mask 直接剔除计算量O(N)O(N²)距离矩阵决策清单两个结构同源、长度一致、要报整体重合度 →RMSD只有预测结构、没有实验参考要评估可信度 →pLDDTAlphaFold 自带见 alphafold/common/confidence.py结构有缺失区域或长度/构象差异大 →lDDT长序列只关心局部折叠质量 →lDDTcutoff 保证它本质是局部指标汇报时把两者都放上注明对齐方式和参与评分的原子集避免单指标误读。踩坑清单对齐方式、缺失原子与长序列几个高频问题按问答收口。Q1RMSD 算出来偏高是不是预测错了不一定。先检查对齐范围——Kabsch 是全局最优旋转末端错 10 Å 也会把全局分数拉高。可以只取匹配残基重算或改用局部 RMSD同时看 lDDT 是否依然很高两把尺子交叉验证。Q2真实结构里有缺失原子PDB 里没解析出来的残基怎么办lDDT 侧把对应位置写进 true_points_mask 即可源码会把它乘进评分掩码自动剔除RMSD 侧则必须手动保证两边选的是同一组残基索引缺的那段两边都跳过去否则对应关系错位分数直接失真。Q3长序列比较cutoff 要调吗15 Å 是局部接触尺度链长不改变逐对比较的语义一般不用动。真正要注意的是全长 RMSD 在长链上会被末端柔化放大建议分域报告而 lDDT 的逐残基分数per_residueTrue能直接指出哪一段接触被破坏比单一总分信息量大得多。Q4lDDT 返回 NaN 或 0哪个是 bug先分清是哪种0残基没有落在 cutoff 内的接触时得分按 1 处理源码用 1e-10 epsilon 兜底见 alphafold/model/lddt.py 第 73 行附近的注释真·0 意味着所有被评分的接触都偏差超过 4 Å那才是真差。NaN 则多半是 mask 或形状没对齐回头查 batch 维和atom_order索引。【免费下载链接】alphafoldOpen source code for AlphaFold 2.项目地址: https://gitcode.com/GitHub_Trending/al/alphafold创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考
返回列表