一、文章介绍理解遗传变异的功能影响是精准医学的核心问题之一。全基因组关联研究GWAS已识别出大量与疾病相关的遗传变异但约90%的疾病相关变异位于非编码区域它们可能通过影响转录或转录后调控发挥作用但其具体分子机制仍不清楚。现有计算预测方法主要分为两类基于进化保守性的方法如PhyloP、PhastCons通过比较多个物种的序列来识别功能受限位点但这些方法通常独立评估每个位点忽略了局部序列模式而基于深度学习的单物种序列模型如DeepSEA、Enformer能够捕捉长程调控信号但不显式利用多物种进化信息。这种割裂导致非编码区域的变异效应预测精度有限——非编码区域缺乏明确的编码约束需要同时依赖局部序列上下文和跨物种进化信号才能准确识别功能相关变异。为弥合这一鸿沟Dongjoon Lim和Mathieu Blanchette在Bioinformatics上发表了题为“Graphylovar: predicting the impact of non-coding variants using a multi-species sequence model”的研究论文提出了一个名为Graphylovar的深度学习方法。Graphylovar的核心创新在于显式地将系统发育树结构整合到深度学习架构中。模型输入为以目标变异为中心的65 bp人类DNA序列以及来自57种其他哺乳动物的直系同源序列和57个计算重建的祖先序列共115个物种/节点。与GPN-MSA等将多序列比对MSA列压缩为单一token而丢弃物种间系统发育结构的方法不同Graphylovar通过图卷积网络GCN沿已知的胎盘哺乳动物系统发育树传播信息使每个物种节点的表征能从其亲缘物种父节点和子节点聚合信息从而保留了进化层次结构。同时Transformer编码器提取每个物种序列内的上下文依赖特征。这种“Transformer 系统发育GCN”的混合架构使模型能从局部序列模式和跨物种进化关系两个维度学习功能约束。Graphylovar采用预测群体水平等位基因频率作为预训练目标而非常见的掩码核苷酸预测。该目标更具生物学信息——在强纯化选择下维持的变异倾向于以较低频率存在于群体中。预训练在TOPMed全基因组测序数据约1.49亿个变异上进行以染色体13-22为完全留出测试集。模型有两个输出头等位基因频率预测头输出A/C/G/T/缺失的五类分布和SNP概率预测头输出该位置为多态性的概率。两者结合为Graphylovar ScoreScoreilog⁡(pspopspr(1−ps))\text{Score}_i \log\left(\frac{p_s p_o}{p_s p_r (1-p_s)}\right)Scorei​log(ps​pr​(1−ps​)ps​po​​)其中po,prp_o,p_rpo​,pr​为模型预测的替代/参考等位基因概率psp_sps​为SNP概率。该分数综合了等位基因特异性约束低概率的替代等位基因更可能有害和位置整体保守性高保守位点的任何改变都更可能有害。在零样本测试约1.49亿留出变异中Graphylovar区分常见变异MAF≥0.01与罕见变异MAF≤0.001的AUROC达0.6246优于PhyloP、PhastCons、CADD、Enformer和GPN-MSA等方法。与CADD的简单集成z-score均值使AUROC提升至0.64420.020P10⁻¹⁵表明Graphylovar捕获了与CADD互补的信息。在13个MPRA基准数据集上微调后的Graphylovar在所有数据集上均取得最高AUROC0.615-0.708优于直接训练的版本和基于Enformer的DNN。物种重要性分析显示人类提供最大的孤立信号贡献而猩猩、绿猴、大猩猩和长臂猿形成第二重要集群非灵长类哺乳动物贡献虽小但广泛存在黑猩猩因与人类序列近乎相同而信息冗余贡献排名仅第11位。二、算法原理介绍Graphylovar的算法设计围绕“Transformer序列编码 → 系统发育GCN传播 → 双头预测与分数合成”三个核心层次展开。一Transformer编码器逐物种序列特征提取。输入为115个物种/节点58个现存物种57个祖先节点的65 bp对齐序列。每个物种的序列经独热编码6字母A/C/G/T/-/N后分别传入正向链和反向互补链两个并行通路。每链先经32维全连接层、批归一化和GELU激活加入正弦位置编码然后通过2层Transformer编码器4个注意力头处理输出经均值池化得到32维表征。正反链的32维表征拼接为64维物种特征向量。Transformer捕获了该物种序列内的短程和长程依赖模式。二GCN与物种注意力系统发育信息传播。64维物种特征向量先通过SESqueeze-and-Excitation块——全局平均池化后经两层全连接和sigmoid产生物种级注意力权重与原始特征逐元素相乘使模型能自动加权重要物种。然后输入2层GCN其邻接矩阵为胎盘哺乳动物系统发育树115个节点边表示进化关系。GCN层更新规则为H(l1)σ(D~−1/2A~D~−1/2H(l)W(l))H^{(l1)} \sigma(\tilde{D}^{-1/2}\tilde{A}\tilde{D}^{-1/2}H^{(l)}W^{(l)})H(l1)σ(D~−1/2A~D~−1/2H(l)W(l))其中A~AI\tilde{A}AIA~AI为带自环的邻接矩阵D~\tilde{D}D~为度矩阵。该操作使每个物种节点从其亲缘物种父母和子女节点聚合进化信息实现跨物种信号的传播与整合——这是Graphylovar与GPN-MSA压缩MSA列的本质区别。三双头预测与Graphylovar Score。GCN输出的115×32维表征经展平3680维后通过共享全连接层128单元ReLU分支为两个预测头1等位基因频率头——64单元softmax输出{A,C,G,T,缺失}上的概率分布2SNP概率头——64单元sigmoid输出该位置为多态性的概率psp_sps​。预训练采用两损失联合优化等位频率头用分类交叉熵SNP头用二元交叉熵。预训练目标为掩码中心人类位置标签来自TOPMed群体等位基因频率。推理时对给定变异参考等位基因rrr替代等位基因ooo模型输出pr,po,psp_r,p_o,p_spr​,po​,ps​。Graphylovar Score定义为log⁡(pspopspr(1−ps))\log\left(\frac{p_s p_o}{p_s p_r (1-p_s)}\right)log(ps​pr​(1−ps​)ps​po​​)。当psp_sps​高可变位点时分数主要由po/prp_o/p_rpo​/pr​驱动低概率替代等位基因给出低有害分数当psp_sps​低高度保守位点时分子pspop_s p_ops​po​变小任何替代都产生低分数正确反映保守位点变异的有害性。三、总结Graphylovar是一个通过显式建模系统发育树来预测非编码变异功能影响的深度学习框架将Transformer提取的逐物种序列特征与GCN沿胎盘哺乳动物系统发育树传播的跨物种进化信号相结合以预测群体等位基因频率为预训练目标并通过双头架构合成兼顾等位基因特异性和位置保守性的Graphylovar Score。其核心创新在于用GCN保留多序列比对中物种间的进化层次结构而非将其压缩为扁平token使模型能从人类近亲到远缘哺乳动物中系统性地学习功能约束。Graphylovar在区分常见/罕见变异AUROC 0.625和MPRA调控效应预测13个数据集均最优上优于现有方法与CADD集成进一步提升了性能且物种重要性分析揭示了符合生物学直觉的进化信息权重分配为非编码变异的功能优先排序提供了兼具预测力和可解释性的新工具。