← 2026-09-15 日报
Nature MethodsVersion of record 更新 重点阅读

SPLENDID:连续祖源交互与稀疏多祖源预测

SPLENDID incorporates continuous genetic ancestry in biobank-scale data to improve polygenic risk prediction across diverse populations

Tony Chen, Haoyu Zhang, Rahul Mazumder, Xihong Lin

论文原文 ↗ PDF ↗ 代码 ↗ Version of record

用基因型×祖源主成分交互刻画连续变化的预测效应,以 group-L0L1 选择变异组,再进行限制变异并集大小的稀疏集成。精读重点是筛选、非凸优化与独立调参的分工,以及预测优势中有多少真正来自交互。

本文目录

从祖源协变量校正转向连续变化的遗传效应

SPLENDID 的目标不是估计每个 SNP 的因果效应,而是用多祖源个体水平数据建立一个不必先给新个体分配离散祖源标签的预测模型。它把“不同祖源是否应有不同预测系数”改写为“预测系数能否随祖源坐标连续变化”,再对这些变系数项做稀疏选择。

这比在常规 PRS 回归中加入几个 PCs 更进一步。PC 主效应只能改变预测的基线;基因型与 PC 的乘积项才能使同一个变异的预测效应随个体祖源坐标变化。不过,不使用离散标签,不等于不使用祖源信息,也不等于在任何未覆盖的人群中都能可靠外推。本文仍以 1000 Genomes 的参考主成分描述遗传背景。来源:主文 Results、Methods 与 Supplementary Methods。

模型:共享主效应与 G×PC 交互

设 GijG_{ij} 为个体 ii 在变异 jj 的基因型,PCikPC_{ik} 为投影到共同参考空间的主成分,ZiZ_i 包含其他协变量。作者的模型可写为

Yi=∑j=1MGijβj0+∑j∈Hpˉ∑k=1K0GijPCikβjk+∑k=1KγkPCik+Zi⊤α+ϵi.Y_i=\sum_{j=1}^{M}G_{ij}\beta_{j0} +\sum_{j\in H_{\bar p}}\sum_{k=1}^{K_0}G_{ij}PC_{ik}\beta_{jk} +\sum_{k=1}^{K}\gamma_kPC_{ik}+Z_i^\top\alpha+\epsilon_i.

本文使用 K=20K=20 个 PCs 做协变量校正,前 K0=5K_0=5 个 PCs 进入交互。对进入异质性候选集的变异,个体对应的预测系数为

Bij=βj0+∑k=15PCikβjk.B_{ij}=\beta_{j0}+\sum_{k=1}^{5}PC_{ik}\beta_{jk}.

不进入候选集的变异只保留共享主效应。这里的 BijB_{ij} 是拟合模型中的个体依赖系数,不是直接测得的个体因果效应。全局祖源 PCs 也不是每个位点的 local ancestry;标签 SNP 的 LD 差异、未建模的因果变异及与祖源相关的环境,都可能影响这种预测系数。

共同 PC 空间是实现的前提。补充方法先在 2,504 名 1000G 个体中进行 LD pruning 和 PCA,再投影训练、调参与验证数据;不能在三个队列分别重新做 PCA 后直接相乘,因为轴的方向、顺序与尺度未必一致。作者选择五个交互 PCs 的依据包括模拟中的祖源效应重构表现,并非对全球遗传多样性的普适维数结论。

异质性筛选不等于最终交互发现

在训练数据中,先逐 SNP 拟合包含主效应、PC 主效应及五个交互项的边际模型,以联合检验

H0:βj1=⋯=βj5=0H_0:\beta_{j1}=\cdots=\beta_{j5}=0

得到异质性 P 值 pjp_j。然后按 pj<pˉp_j<\bar p 决定是否为该 SNP 开放交互项。补充方法使用

pˉ∈{0, 5×10−8, 5×10−6},\bar p\in\{0,\ 5\times10^{-8},\ 5\times10^{-6}\},

其中零阈值提供不含交互的候选模型。这一步同时缩小计算规模,并让后续预测选择在“全共享”和“部分异质”之间调节。

需要区分两个任务:筛选 P 值决定进入哪类预测模型,最终是否保留及如何加权还由惩罚回归和调参数据决定。因此,最终选中的 G×PC 项不能直接列为经过全流程错误率控制的生物学交互。边际筛选也可能遗漏在联合条件分析后才显现的异质性;这是计算筛选与模型完整性之间的取舍。

Group-L0L1:选择变异组,而不是逐个交互项做 Lasso

对每个 SNP 定义一组系数。未通过筛选时为标量 βj0\beta_{j0},通过时为

βj=(βj0,βj1,…,βj5)⊤.\boldsymbol\beta_j=(\beta_{j0},\beta_{j1},\ldots,\beta_{j5})^\top.

作者最小化残差平方和与分组惩罚之和:

L=∑i(Yi−Y^i)2+∑j{λ01(∥βj∥2>0)+λ1∥βj∥2}.\mathcal L=\sum_i(Y_i-\widehat Y_i)^2+ \sum_j\left\{\lambda_0\mathbf1(\|\boldsymbol\beta_j\|_2>0) +\lambda_1\|\boldsymbol\beta_j\|_2\right\}.

L0L_0 项计数被保留的变异组,L1L_1 项则是组间求和、组内取 L2L_2 范数的 group Lasso,而不是对所有系数分别取绝对值。这样可以一方面限制参与预测的 SNP 数量,另一方面减小选中组内系数的方差。

该写法把主效应与交互放入同一组,但不额外施加“只要存在交互,主效应系数就必须非零”的强遗传层级约束。稀疏选择也不产生 SuSiE 意义下的 PIP 或 credible set。LD 中选出一个高预测价值的代表 SNP,与正确确定分子因果变异,是不同的成功标准。

算法:组坐标更新、硬阈值与 active set

补充方法给出了具体计算步骤。令 AjA_j 为 SNP jj 的设计矩阵块,包含 GjG_j 以及获准进入的 GjPCkG_jPC_k;令 LjL_j 为 Aj⊤AjA_j^\top A_j 的最大特征值,rr 为当前残差。在作者未除以样本量的平方损失约定下,先计算

uj=βj+Aj⊤rLj.u_j=\boldsymbol\beta_j+\frac{A_j^\top r}{L_j}.

再作组更新:

βjnew=uj(1−λ12Lj∥uj∥2)1{∥uj∥2>λ0Lj+λ12Lj}.\boldsymbol\beta_j^{\rm new}= u_j\left(1-\frac{\lambda_1}{2L_j\|u_j\|_2}\right) \mathbf1\left\{\|u_j\|_2> \sqrt{\frac{\lambda_0}{L_j}}+\frac{\lambda_1}{2L_j}\right\}.

式中的 uju_j 即上式定义的更新向量;零向量不进入保留分支。这个操作把组内收缩与整组删除结合起来。算法循环更新活跃组,达到收敛条件后再检查未激活组是否需要纳入。它利用分组二次近似与 active-set 策略获得高质量解,不是对非凸最优解的全局保证。实际计算还借助 bigstatsr/bigsnpr 的内存映射数据结构,避免简单展开全部大矩阵所需的内存。

补充方法的参数路径包括五个对数尺度 λ1\lambda_1 值(10−210^{-2} 到 10−410^{-4})及每个路径最多 50 个 λ0\lambda_0 值;外部调参误差不再下降或支持集达到上限时可提前停止。因此,预测模型与运行时间同时依赖搜索范围、停止规则和数据布局,而非只依赖优化公式。

稀疏集成与没有外部调参集时的 CMSA

主流程分为独立的训练、调参和验证数据。训练阶段产生多组 (λ0,λ1,pˉ)(\lambda_0,\lambda_1,\bar p) 的 PRS;调参阶段先选出性能最佳单模型,其变异集合记作 SS。再按候选模型表现依次纳入其他 PRS,但限制所用变异并集小于 c∣S∣c|S|,本文 c=2c=2;最后通过交叉验证 ridge regression 组合这些分数。

所以调参集并非重复训练一次 GWAS,而是在选择正则化强度、是否需要交互以及模型组合权重。即使原始遗传系数不再直接用调参表型拟合,模型选择也已使用其表型,必须留下完全独立的最终测试集。

作者还给出 CMSA:在 AoU 内部五折训练,每一折用四折拟合、一折选择最佳模型,再平均各折系数。主文 Extended Data Fig. 5 显示,它与使用 UKB 外部调参的 SPLENDID 表现接近。这是对调参迁移稳健性的有用检查,但不是“任何分布偏移都不会过拟合”的证明。

外部 GWAS 有两种不同用途:以 Cochran’s Q 替换 G×PC 筛选,以及利用外部效应对惩罚作自适应加权。外部 GWAS 排除 UKB 有助于避免训练—测试重叠。补充公式的自适应权重写法仍有符号处理不清之处,本文不将其补写成未经实现核验的完整公式。

实证设计:队列总规模不能替代性状样本量

层次原文设计及其含义
训练AoU 可用无关个体 224,364 人;使用阵列 v7.1,不是将所有 WGS 变异输入模型
调参与验证UKB 340,140 人随机分半,约各 17 万;最终非欧洲祖源验证约 1.2 万
变异AoU、UKB 插补数据和 1000G 共有的 1,172,070 个变异;至少一个 AoU 预测祖源组 MAF>1%
表型身高、四项血脂及四项血细胞性状;提取最近测量,使用他汀者 LDL-C/0.7、TC/0.8
性状有效样本补充表 5 中 AoU LDL-C 为 89,211、HDL-C 为 94,165、log-TG 为 95,451、身高为 222,399;不能统称各性状均有 22 万训练样本
评估从 PRS 和表型中去除年龄、性别及 20 PCs 后比较 partial R²;另给出 MSE 类指标及 10,000 次 bootstrap

预测时,SPLENDID 的 PRS 包含遗传主效应和 G×PC 项,不包含训练队列中估计的年龄、性别和 PC 主效应。模型在拟合时调整这些协变量,与把协变量本身计入 PRS,是两回事。

**公式核对提醒:**主文 Methods 的显示公式把 Rpartial2R^2_{\rm partial} 写成了未平方的 cor,而正文、图表与补充表均以 partial R² 表述。本笔记沿用作者对性能指标的名称,同时保留这一印刷表达不一致;不能把该等式原样当作相关系数与决定系数等价的数学定义。

性能收益有多大,又来自哪里

九性状非欧洲祖源验证中,作者报告平均 partial R² 相对 PRS-CSx、CT-SLEB、PROSPER 和 iPGS+refit 分别提高 65.8%、42.6%、29.4% 和 34.7%。这里均是相对提高。补充表 7A 的两个例子使尺度更清楚:

非欧洲祖源验证SPLENDID partial R²(95% CI)PRS-CSx partial R²(95% CI)
LDL-C0.068(0.059–0.077)0.041(0.035–0.049)
HDL-C0.098(0.087–0.110)0.081(0.071–0.092)

主文 Figure 3 应结合这个绝对尺度阅读。更大的相对比例,可能对应基线较低时有限的绝对增加。

消融尤其重要:补充表 11 中 LDL-C 的无交互模型仍约为 0.068,加入交互的相对差异为 0.1%,置信区间跨零。相反,主文报告 MCH 和 NEUT 的交互纳入分别带来约 5% 和 20% 的收益。因此,SPLENDID 的优势由合并个体数据、组稀疏正则化、集成以及在部分性状中有效的交互共同组成,不能只归因于“祖源连续化”。

它也不是所有比较的赢家。欧洲祖源中有七个性状略逊于 CT-SLEB/PROSPER;减少欧洲调参样本后这种格局改变,说明方法比较还依赖调参数据分配。模拟中相对 summary-statistics 方法的优势通常比相对 iPGS+refit 更明显,部分设置后者表现更好。

补充表 9 中,SPLENDID 平均保留 1,803 个 SNP,iPGS+refit 为 3,959。LDL-C 保留 173 个主效应变异,其中 5 个带交互;身高则达到 9,999 个主效应变异,接近实际搜索所设的 10,000 上限。因此,“选中数量接近真实支持集”主要是特定模拟结果,不能据真实身高的稀疏模型推断其真实多基因性只有一万个变异。PC 本身也需要额外的大量基因型标记。

生物学解释与计算基准的边界

作者讨论了 Duffy-null 变异在祖源间 MAF 不同,以及 MCH 中模型预测效应随 PC 变化的例子。这些解释有遗传学依据,但主文写作存在一个应保留的差异:MCH 例子正文为 rs1211375,Figure 4 图注为 rs1213375。在没有进一步核对位点文件前,不能自行确定哪一个是正确标识。

Figure 4 的个体效应分布来自拟合的 BijB_{ij},其祖源组间差异检验是对该模型输出的比较,不是独立实验测量证实了生物学互作。标签变异在不同 LD 背景中得到不同最佳线性系数,也可能改善预测而不意味着因果效应发生同等程度改变。

计算方面,补充表 3 的模拟基准使用 10 万训练样本与约 95 万变异:G×PC 扫描和初始统计量分别约 58–70 和 80–100 分钟,且是多性状共享/并行设置;后续惩罚回归每个性状约 73–185 分钟。不同阶段使用 10–20 个核心,内存峰值约 30 GB。不能把“约五小时”改写为任意 biobank、单核、全流程固定耗时,也不能与 GPU 内核加速倍数直接比较。

比较实现也不是全部未经修改的默认版本:CT-SLEB 的部分阈值被删减,PROSPER 去掉了数值表现失败的 OLS learner,个体水平模型设置了支持集搜索上限。作者说明了这些决定;复现时应保持记录,而不是只安装包名就期待相同结果。

对跨祖源与细胞状态 QTL 方法的启示

以下是由本文引出的研究设想,不是本文已经验证的 QTL 方法。

值得迁移的是共享效应加低维连续偏离的结构,而不是先把细胞状态切成多个箱子再比较显著性。对重复测量的 donor–cell 数据,可考虑让 βj\beta_j 随 donor 祖源坐标与细胞状态坐标变化,同时单独处理 donor 内相关。但若直接同时开放祖源、状态及三阶交互,参数量与不稳定性会迅速增加;需要先比较主效应、单一交互、共享低秩交互等有明确层级的模型。

评估应把两个目标拆开。预测研究需要按 donor、队列留出,检查外部均方误差、校准和连续祖源空间中的最弱表现区域。推断研究则需要重新建立交互检验与多重比较控制,不能拿预测调参后的非零系数作为发现列表。若进一步进入 fine-mapping,功能先验、LD 和效应异质性的建模也应独立校准。

最有信息量的基线不是只有“传统 PRS 对 SPLENDID”,而是相同数据、相同搜索预算下的共享效应 group-L0L1、含交互 group-L0L1、Lasso、不同集成策略与不同 PC 表示。这才能定位收益究竟来自连续背景建模还是更好的正则化与信息利用。

阅读评价与资料

**评价:**这是一篇值得精读的预测方法论文,贡献是把可扩展的分组非凸选择、连续祖源交互和稀疏集成组合成统一流程。它对多祖源个体数据的利用有实证支持,消融与 CMSA 增强了说服力;对因果变异解释、全球祖源外推以及临床绝对风险的证据则明显有限。最值得借鉴的是模型结构和验证拆分,而不是把稀疏输出直接解释为生物学发现。

已核读主文与 Methods、Extended Data 1–5、补充方法的 PC/筛选/优化/CMSA 部分,以及补充表 3、5、7A、9、11 等关键表。其余补充分析按原文报告描述,未在个体数据上重跑。

期刊图像为出版商专有许可;本页提供对应原图链接,不重新分发该文图像或全文 PDF。