MethylGPT: DNA甲基化组基础模型.

paper：MethylGPT: a foundation model for the DNA methylome

0. TL; DR

MethylGPT 是一个基于 transformer 架构的DNA甲基化基础模型，在一个包含超过15万个人类甲基化谱的庞大数据集上进行了预训练，该数据集涵盖了多种组织类型、49,156个关键 CpG 位点和76亿个训练 token。通过无监督学习，MethylGPT 能够学习到 CpG 位点的生物学意义表示，捕捉局部基因组上下文和高阶染色体特征。

MethylGPT 在预测缺失的甲基化值方面表现出色（Pearson 相关系数 $R=0.929$）。在下游任务中，即使有多达70%的数据缺失，模型性能依然稳定。在跨多种组织的年龄预测任务中，其准确性优于现有方法。在对60种主要疾病进行风险预测时表现稳健，并能够系统地评估干预措施对疾病风险的影响。

1. 背景介绍

DNA 甲基化是一种关键的表观遗传修饰，通过在 CpG 二核苷酸的胞嘧啶残基上添加甲基基团来调控基因表达。它不仅在细胞分化、发育和维持基因组稳定性等基础生物学过程中扮演核心角色，还因其稳定性、对环境因素的动态响应以及在临床样本中的可及性，成为一种理想的生物标志物。

近年来，全基因组甲基化谱分析已经在多种病理状态（特别是癌症检测和心血管风险评估）中揭示了独特的特征。同时，与年龄相关的甲基化模式催生了高精度的表观遗传时钟，如 DunedinPACE 和 GrimAge，这些时钟在预测年龄相关疾病和寿命方面甚至优于传统的临床指标。

然而，将基于甲基化的诊断方法应用于临床仍面临诸多分析挑战。目前的计算方法主要依赖于线性模型（如 ElasticNet）和简单的统计方法。这些模型的根本局限性在于：

无法捕捉非线性关系：它们难以模拟甲基化数据中复杂的、非线性的相互作用。
忽略位点间的依赖性：它们通常假设 CpG 位点之间相互独立，忽略了构成甲基化模式的调控网络和高阶相互作用。
缺乏上下文感知：相同的甲基化模式在不同的细胞和组织背景下可能具有不同的生物学意义，这是线性模型无法捕捉的。
对技术噪音敏感：线性模型难以处理批次效应和缺失数据等引入的非线性变异。

与此同时，人工智能领域的进展，特别是 transformer 架构和基础模型，已经彻底改变了复杂生物序列的分析。从蛋白质组学（ESM-2、AlphaFold3）到基因组学（Enformer、Evo），再到单细胞领域（Geneformer、scGPT），基础模型在学习跨任务泛化的通用生物学模式方面展现了非凡的能力。然而，DNA 甲基化分析领域长期缺乏这样一个统一的框架。

基于此，作者引入了 MethylGPT，一个专为 DNA 甲基化设计的基于 transformer 的基础模型，旨在提供一个统一的分析框架，既能保留生物学上下文，又能适应多样化的下游任务。

2. MethylGPT 方法

作者构建 MethylGPT 的核心思想是利用 transformer 架构学习 DNA 甲基化谱中的复杂依赖关系。整个流程包括数据收集与预处理、模型架构设计、预训练和下游任务微调。

2.1 数据收集与预处理

作者从 EWAS Data Hub 和 Clockbase 中收集了226,555个人类 DNA 甲基化谱。经过严格的质量控制（quality control, QC），包括过滤掉缺失率超过40%的低质量样本，并进行去重（deduplication），最终保留了154,063个高质量样本用于预训练。这些样本覆盖了超过20种不同的组织类型。

为了确保生物学相关性和跨平台兼容性（Illumina 27k, 450k, EPIC 芯片），作者筛选了49,156个生理学上重要的 CpG 位点。筛选标准为：该 CpG 位点在 EWAS 目录中与超过5个性状相关联，且在超过95%的预训练数据集中出现。

最终，数据被处理成一个矩阵 $X \in R^{N \times M}$，其中 $N$ 是样本数量，$M$ 是 CpG 位点数量（49,156）。每个元素 $x_{i,j}$ 代表样本 $i$ 中 CpG 位点 $j$ 的甲基化值（beta-value）。

2.2 模型架构

MethylGPT 由三个主要部分组成：嵌入模块、transformer 模块和任务特异性头（task-specific heads）。

2.2.1 CpG 位点令牌化（Tokenization）

对于49,156个 CpG 位点中的每一个，作者都分配了一个唯一的整数标识符 $id \in {0, 1, …, M-1}$。

2.2.2 嵌入层（Embedding Layers）

对于每个 CpG 位点，模型不仅编码其身份，还编码其甲基化状态。具体来说，嵌入模块将 CpG 位点的标识符和其对应的甲基化值分别投射到两个独立的嵌入向量中，然后通过逐元素求和（element-wise sum）将它们合并。

对于样本 $i$ 的最终嵌入 $E_{\text{final}}$ 定义为：

\[E_{\text{final}}^{(i)} = E_{\text{CpG ID}} + E_{\text{Methylation Value}}^{(i)}\]

其中，$E_{\text{CpG_ID}}$ 是所有样本共享的 CpG 位点身份嵌入，$E_{\text{Methylation_Value}}^{(i)}$ 是样本 $i$ 特有的甲基化值嵌入。嵌入维度设置为64。此外，在每个序列的开头，作者添加了一个特殊的 [CLS] token，用于学习样本级别的全局表示。

2.2.3 Transformer 模块

作者采用了包含6个 transformer 模块的标准自注意力（self-attention）编码器。每个模块包含一个多头自注意力层（4个头）和一个标准的前馈网络（MLP）。每个子层之后都应用了层归一化和残差连接。

transformer 的处理过程可以表示为：

\[\begin{align*} h'^{(l)} &= \text{MSA}(h^{(l)}) + h^{(l)} \\ h^{(l+1)} &= \text{MLP}(h'^{(l)}) + h'^{(l)} \end{align*}\]

为了高效地处理长达数万个 CpG 位点的序列，作者利用了 FlashAttention 的实现，它在加速训练和推理的同时，最大限度地减少了内存占用。

2.2.4 任务特异性头

在预训练阶段，一个简单的线性层将 transformer 的输出嵌入投射回单个值，以预测甲基化水平。在微调阶段，作者使用更复杂的网络（如 ResNet1D）来处理完整的输出嵌入，以预测生物年龄或疾病风险。

2.3 模型预训练

模型使用 AdamW 优化器进行训练，学习率为0.001，批大小为16，在 NVIDIA A100 GPU 上训练了10个周期。

预训练阶段采用了两个互补的目标函数：

甲基化值预测（Methylation Value Prediction, MVP）

似于掩码语言模型（Masked Language Modeling, MLM）。作者随机掩盖（masked）30%的 CpG 位点（即不输入它们的甲基化值嵌入），然后训练模型来预测这些被掩盖位点的原始甲基化值。目标是最小化均方误差（Mean Squared Error, MSE）。

\[\mathcal{L}_{\text{MVP}} = \frac{1}{|\mathcal{M}|} \sum_{j \in \mathcal{M}} (\hat{y}_j^{(i)} - y_j^{(i)})^2\]

其中 $\mathcal{M}$ 是被掩盖位点的集合，$\hat{y}$ 是预测值，$y$ 是真实值。

谱重构（Profile Reconstruction）

利用 [CLS] token 的输出嵌入来重构完整的甲基化谱。模型将 [CLS] token 的输出嵌入反馈回输入，同时掩盖所有其他 CpG token。目标同样是最小化 MSE。

\[\mathcal{L}_{\text{Reconstruction}} = \frac{1}{|\mathcal{A}|} \sum_{j \in \mathcal{A}} (\hat{y}_j^{(i)} - y_j^{(i)})^2\]

其中 $\mathcal{A}$ 是所有位点的集合。

3. 实验分析

作者进行了一系列详尽的实验来验证 MethylGPT 的性能、生物学相关性以及在下游任务中的应用潜力。

3.1 模型开发与验证

图 b: 训练曲线显示，模型在训练过程中收敛迅速且没有明显的过拟合。在第10个周期时，模型在测试集上达到了最佳的均方误差，为0.014。
图 c, d, e, f: 在缺失值插补（imputing）任务中，MethylGPT 表现出色。预测值与真实值之间的 Pearson 相关系数高达 $R=0.929$，平均绝对误差为0.074。残差图和分段误差分析显示，模型在不同甲基化水平（从0.0到1.0）上都能保持稳健的预测精度。

3.2 MethylGPT学习到具有生物学意义的CpG表示

作者使用 UMAP 对49,000个 CpG 位点的上下文嵌入进行降维可视化，以探究模型是否在无监督的情况下学到了生物学特征。

图 b: UMAP 可视化显示，CpG 位点根据其与 CpG 岛（island）的关系（岛、岸shore、架shelf和其他区域）形成了清晰的聚类。
图 c: 嵌入空间中，增强子（enhancer）区域的 CpG 位点也显示出独特的聚类模式。
图 d: 嵌入清楚地将性染色体（sex chromosomes）与常染色体（autosomes）分离开来。

这些结果表明，MethylGPT 成功地捕捉了从局部序列上下文到高阶染色体特征的多层次调控信息，而这些信息是传统独立位点分析方法无法获得的。

3.3 MethylGPT学习到组织和性别特异性甲基化模式

接下来，作者分析了样本级别的嵌入，以评估模型是否能捕捉到有意义的生物学分组。

图 a, b, c: 使用 MethylGPT 生成的样本嵌入进行 UMAP 可视化后，样本根据组织类型（如全血、大脑、肝脏）和性别形成了非常清晰的聚类。同时，批次效应（batch effects）在嵌入空间中不显著。
图 d, e, f: 作为对比，作者直接使用原始甲基化数据进行 UMAP 分析。结果显示，原始数据的聚类效果远不如 MethylGPT 的嵌入：组织边界模糊，批次效应更强，性别分离也较弱。

这有力地证明了 MethylGPT 的上下文嵌入方法能够放大生物学相关信号，同时抑制技术噪音。

3.4 MethylGPT在年龄预测中表现出色

作者在一个包含11,453个多组织样本的数据集上对 MethylGPT 进行了年龄预测的微调。

图 b, c, d: 在微调之前，预训练的嵌入已经显示出与年龄相关的组织结构。经过微调后，这种年龄依赖的聚类变得更加明显，同时组织特异性模式也得以保留。
图 e: MethylGPT 的年龄预测性能优于所有基准方法，包括 ElasticNet、MLP (AltumAge) 和 Horvath 的皮肤和血液时钟。在验证集上，MethylGPT 实现了4.45年的中位绝对误差（Median Absolute Error, MedAE）。
图 f: 在鲁棒性分析中，MethylGPT 表现出惊人的韧性。即使在高达70%的数据缺失情况下，其性能依然保持稳定，显著优于 ElasticNet 和 MLP。
图 g, h: 作者将模型应用于 iPSC 重编程过程中的甲基化数据。MethylGPT 的嵌入清晰地揭示了一个返老还童（rejuvenation）的轨迹。其预测的表观遗传年龄在重编程过程中急剧下降，与 GrimAge 和 Horvath 时钟的观察结果一致，这为 iPSC 重编程是一种真正的衰老逆转方法提供了独立证据。

3.5 年龄特异性注意力模式揭示了不同的甲基化特征

为了深入理解 MethylGPT 是如何处理年龄相关信息的，作者分析了模型的注意力权重。

图 a, b: 模型为年轻（<20岁）和年老（>60岁）样本学习到了截然不同的 CpG 位点交互模式。注意力权重集中在少数 CpG 位点上，表明这些稀疏的位点可能对年龄特异性甲基化至关重要。
图 c, d: 通过火山图识别出在年轻和年老组之间具有显著差异注意力的 CpG 位点。热图分析显示：
- 在年轻样本中，高注意力位点主要与非年龄相关的表型（如性别、自身免疫病）相关。
- 在年老样本中，高注意力位点则与衰老及衰老相关性状（如 BMI、甲状腺病变）密切相关。
图 e: 功能富集分析进一步证实了这些观察。
- 年轻组的高注意力位点富集在发育过程中，如对生长因子刺激的细胞反应。
- 年老组的高注意力位点则富集在氧化应激和氨基酸代谢等与衰老相关的通路中。

这些结果表明，MethylGPT 的注意力机制能够捕捉到具有生物学意义的、随年龄变化的甲基化调控模式。

3.6 疾病风险预测与干预分析

最后，作者在 Generation Scotland 队列（n=18,859）上评估了 MethylGPT 的临床应用潜力。

图 a, b, c: 作者微调模型以预测60种主要疾病（分为8大类）在10年内的风险。模型在验证集和测试集上分别取得了0.74和0.72的总体 AUC，展示了其稳健的预测能力。
图 d: 作者利用这个疾病预测框架，系统评估了8种不同干预措施对预测疾病风险的影响。热图显示了干预措施特异性的效果：
- 戒烟（Smoking cessation） 对降低10年死亡率和心血管疾病风险具有最强的保护作用。
- 高强度训练（High-intensity training） 对呼吸系统、神经系统和自身免疫性疾病有显著益处。
- 依维莫司（Everolimus） 治疗显示会增加自身免疫性疾病的风险，这与之前关于长期使用免疫抑制剂的研究结果一致。

这些发现证明了 MethylGPT 在预测干预特异性健康结果和优化个性化干预策略方面的巨大潜力。