CpGPT: DNA甲基化基础模型.

paper：CpGPT: a Foundation Model for DNA Methylation

0. TL; DR

Cytosine-phosphate-Guanine Pretrained Transformer (CpGPT) 是一个DNA 甲基化基础模型，在超过1,500个DNA甲基化数据集、涵盖超10万个来自不同组织与条件的样本上进行预训练，利用改进的Transformer架构学习甲基化模式的综合表征，能够从有限的输入数据中插补和重建全基因组甲基化图谱。

通过捕获序列、位置和表观遗传上下文，CpGPT在针对衰老相关任务（如时序年龄预测、死亡风险和发病率评估）进行微调时，性能优于专业模型。该模型能高度适应不同的甲基化平台和组织类型。此外，通过分析样本特异性注意力权重，能够识别对个体预测最具影响力的CpG位点。

1. 背景介绍

DNA 甲基化，特别是在 CpG 二核苷酸上胞嘧啶的甲基化，是基因表达表观遗传调控的核心机制，深刻影响细胞分化、组织发育、衰老以及众多疾病（如癌症、神经退行性疾病）的进程。随着高通量甲基化芯片（如 Illumina 450K, EPIC）和测序技术的发展，海量甲基化组数据得以产生。

传统分析方法，如基于线性回归的表观遗传时钟（例如 Horvath时钟、Hannum时钟、GrimAge），虽然成功地将甲基化模式与生物学年龄、死亡率等表型关联起来，但其模型假设相对简单，难以捕捉CpG位点间复杂的非线性相互作用和长程依赖关系。此外，不同甲基化平台之间的数据整合、缺失值的插补、以及模型本身的可解释性，一直是该领域的挑战。

近年来，Transformer架构及其在自然语言处理（NLP）和计算生物学（如 scGPT、DNA语言模型）中的成功应用，为处理序列化和结构化的生物医学数据提供了新思路。将Transformer的强大表征学习能力应用于甲基化组数据，有望学习到更通用、更复杂的表观遗传规律，从而成为一个强大的基础模型，服务于下游多种预测和解释任务。

为此，作者开发了 CpGPT，一个专门为DNA甲基化设计的基础模型，旨在统一处理上述挑战，并为衰老和疾病研究提供新的分析工具。

2. CpGPT 方法

CpGPT 的核心设计思想是将每个CpG位点视为一个词元，其序列背景、基因组位置和甲基化状态（β值）作为输入特征，利用Transformer架构学习其上下文感知的嵌入表示。

2.1 数据与预处理

作者构建了名为 CpGCorpus 的大规模数据集，整合了来自 GEO 数据库的 1,502 项研究、共 106,795 个人类样本。数据涵盖多种 Illumina 平台（27k, 450k, EPIC, EPIC+, EPICv2）和多样化的组织、疾病状态。

预处理步骤包括：

使用 SeSAMe 流程处理原始 IDAT 文件。
对已处理数据，进行质量控制，过滤低质量样本和探针。
探针协调：基于基因组坐标将不同平台的探针映射到统一的CpG位点集。

数据集被划分为训练集（100,965样本）、验证集（489样本）和测试集（5,341样本）。

2.2 模型架构

CpGPT 的输入是每个CpG位点的三元组信息：序列嵌入 $E$、甲基化$β$值 $\beta$、位置信息 $p$。

作者训练了两个版本的模型：一个大型模型（65M参数）用于零样本任务，一个小型模型（2.5M参数）用于微调。

序列编码

对于每个CpG位点，提取以其为中心、长度为 $L=2001$ 的核苷酸序列(-1000bp~1000bp)。使用预训练的DNA语言模型（如 Nucleotide Transformer）得到序列嵌入 $E \in \mathbb{R}^{b \times l \times d_{dna}}$，其中 $b$ 为批次大小，$l$ 为序列长度，$d_{dna}$ 为DNA嵌入维度。随后通过一个可训练的多层感知机 $f_{\text{seq}}$ 将其投影到模型内部嵌入空间：

\[S = f_{\text{seq}}(E)\]

其中 $S \in \mathbb{R}^{b \times l \times d_{\text{emb}}}$，$d_{\text{emb}}$ 为模型嵌入维度。

双重位置编码

为捕获多尺度基因组位置信息，采用了双重编码策略：

绝对位置编码：适应基因组距离的正弦编码。令 $p_{\text{max}} \approx 3 \times 10^8$（最大人类染色体长度），生成位置编码 $PE(p) \in \mathbb{R}^{b \times l \times d_{\text{emb}}}$：

\[t = \frac{p}{p_{\text{max}}}, \quad w = 2\pi \frac{p}{p_{\text{max}}}\]

相对位置编码：应用旋转位置嵌入（RoPE）来高效捕获局部位置关系：

\[L = \text{RoPE}(S + PE(p))\]

甲基化状态编码

甲基化$β$值通过一个可训练的MLP $f_\beta$ 进行嵌入：

\[B = f_\beta(\beta)\]

最终每个CpG位点的嵌入 $C$ 由位点嵌入和$β$值嵌入组合并归一化得到：

\[C = \frac{1}{\sqrt{2}} (L + B)\]

Transformer++

CpGPT 使用 Transformer++ 架构作为核心，其改进包括SwiGLU激活函数、RMSNorm、前置层归一化、无偏置线性层。

模型包含 $N$ 层 Transformer++ 块。在输入序列前添加一个可学习的分类词元 $\mathbf{c}_{\text{cls}}$，用于聚合序列信息，得到最终样本嵌入 $\mathbf{h}_{\text{cls}}$。

解码器

CpGPT设计了多个专用解码器：

$β$值解码器：预测CpG位点的甲基化$β$值。 $\hat{\beta} = \sigma(\langle g_\beta(S), \mathbf{h}_{\text{cls}} \rangle)$
不确定性解码器：估计预测的不确定性。 $\hat{\sigma}_\beta = \sigma(\langle g_{\text{unc}}(S), \mathbf{h}_{\text{cls}} \rangle)$
条件解码器：用于下游任务的微调。 $\hat{y} = \langle g_{\text{cond}}(T), \mathbf{h}_{\text{cls}} \rangle$

2.3 预训练与损失函数

对每个训练批次，将$β$值按照$50\%$随机分成输入集 $\beta_{\text{input}}$ 和目标集 $\beta_{\text{target}}$，鼓励模型重建缺失数据。

预训练采用多任务学习。损失项包括：

$β$值损失：同时采用平均绝对误差（MAE）和Wasserstein距离损失，MAE 损失引入了不确定度预测：

\[\mathcal{L}_{\text{mae}} = \frac{1}{|\mathbf{M}_{\text{valid}}|} \sum_{i,j \in \mathbf{M}_{\text{valid}}} |\hat{\beta}_{ij} - X_{ij}| \\ \mathcal{L}_{mae\_unc} = \frac{1}{|\mathbf{M}_{valid}|}\sum_{i,j \in \mathbf{M}_{valid}} ||\hat{\beta}_{ij} - X_{ij}| - \hat{\sigma}_{\beta,ij}| \\ \mathcal{L}_{\text{wd}} = \mathbb{E}[|\text{sort}(\hat{\beta}) - \text{sort}(X)|]\]

样本嵌入损失：KL散度损失约束嵌入空间为正态分布；对比损失拉近相似样本、推开不相似样本。

\[\mathcal{L}_{kld\_sample} = \text{KL}(\mathcal{N}(\mathbf{\mu}_h, \sigma_h^2) \| \mathcal{N}(0, \mathbf{I})) \\ \mathcal{L}_{contrastive} = \mathbb{E}[1 - (\max(0, \cos(\mathbf{h}_{sample,i}, \mathbf{h}_{sample,j}) - \tau))^2]\]

其中 $\tau$ 为阈值参数，$\mu_h, \sigma_h^2$ 为 $\mathbf{h}_{sample}$ 的批次均值和方差。

条件预测损失：根据下游任务类型（回归、分类、生存分析）选择特定损失函数，例如用于生存分析的Cox比例风险损失：

\[\mathcal{L}_{CPH} = -\frac{1}{\sum_{i=1}^b e_i} \sum_{i:e_i=1} \left(\hat{y}_i - \log\left(\sum_{j:t_j \geq t_i} \exp(\hat{y}_j)\right)\right)\]

其中 $\hat{y}_i$ 为预测风险评分，$t_i$ 为观察时间，$e_i$ 为事件指示符。

3. 实验分析

作者对 CpGPT 进行了系统性的评估，涵盖了表征学习、零样本任务、可解释性以及在下游衰老和疾病预测任务上的性能。

3.1 CpGPT学习到有意义的嵌入

为了验证模型是否无监督地学习到了有意义的表征，作者对CpG位点嵌入和样本嵌入进行了 UMAP 可视化分析。

图a-b：CpG位点嵌入的UMAP图显示，模型能根据CpG岛注释和ChromHMM染色质状态注释自然地将位点聚类，表明模型从数据中自发捕获了功能基因组注释信息。
图c-d：样本嵌入的UMAP图在 AltumAge 数据集和 RRBS Atlas 数据集中清晰地按组织类型分离，例如血液、大脑、肝脏、皮肤等组织形成不同的簇。
图e-f：零样本参考映射实验。将 Hannum 数据集（GSE40279）的样本映射到 AltumAge 数据的组织注释上。结果显示，尽管 Hannum 数据主要为全血样本，CpGPT 能准确地将大部分样本映射到正确的血液相关亚型（如全血、白细胞），证明了其样本嵌入具有强大的生物学信息表征能力和零样本泛化能力。

CpGPT 的嵌入空间无需显式监督，即可编码丰富的基因组功能语境和样本生物学状态信息。

3.2 有效的零样本插补与平台转换

CpGPT 被设计用于从部分数据重建全甲基化组。

图a-b：在 Hannum 数据集上模拟从 MSA 芯片数据重建 450k 芯片探针。对于MSA芯片缺失的450k探针，CpGPT 插补的$β$值与真实值高度相关（Pearson 相关系数中位数达 0.991，Spearman 相关系数中位数 0.964，MAE 中位数 0.032）。
图c-d：使用插补后的450k数据预测多种表观遗传时钟（如 AltumAge、GrimAge、PACE 等），其预测值与使用真实450k数据得到的结果高度一致，且性能显著优于直接使用原始MSA有限探针的预测。热图显示，对于几乎所有测试的生物标志物，CpGPT 重建数组的标准化MAE都低于原始MSA数据。
图e：将模型在多种哺乳动物数据上微调后，对 12 种未见过的物种进行Horvath哺乳动物芯片一半位点的插补，结果显示，与人类亲缘关系更近的物种（如大猩猩、猕猴）插补性能最佳（Pearson 相关接近 0.9），证明了模型跨物种的泛化能力。

CpGPT 能够高精度地插补缺失的甲基化值，实现不同甲基化平台间的有效转换，解决了表观遗传研究中数据整合的关键难题。

3.3 样本特异性的甲基化组解读

Transformer 的注意力机制为模型提供了内在的可解释性。作者通过分析最后一个Transformer块的平均注意力权重，为每个样本构建了CpG位点-位点重要性矩阵。

作者选取了 AltumAge 数据集中心脏组织和大脑组织的样本进行分析，仅保留注意力权重最高的前 0.01% 的连接构建网络图。

心脏样本：连接最紧密的CpG位点位于 11 号染色体的 KCNQ1 基因内。KCNQ1 编码心脏动作电位复极化关键的电压门控钾通道，其突变可导致长QT综合征。模型自动识别出该基因在心脏组织中的核心地位。
大脑样本：连接最紧密的CpG位点位于 5 号染色体的 GDNF 基因内。GDNF 是一种重要的胶质细胞源性神经营养因子，对神经元存活、发育和可塑性至关重要。

CpGPT 的注意力机制能够提供样本特异性的解读，动态地突出与特定组织生物学功能最相关的基因和CpG位点，这超越了传统线性模型赋予固定权重的静态解读方式。

3.4 死亡率与发病率预测能力

作者在多个独立队列上微调 CpGPT 用于死亡率和发病率预测，并采用Cox比例风险模型进行评估。

死亡率预测：

图A-C：在三个训练队列和一个测试队列中，CpGPT 预测死亡率的C-index在训练队列中范围为 0.68 到 0.82，在测试队列中达到 0.82。AUC值在训练队列中为 0.70 到 0.87，在测试队列中高达 0.90。z-score也显示了预测因子与死亡率的强关联。
图D-F：根据年龄残差化的CpGPT风险评分将个体分为加速衰老（最高四分位）和减速衰老（最低四分位）组。Kaplan-Meier生存曲线显示，在所有队列中，两组间的生存概率存在极显著差异（$p < 0.0001$）。在测试队列中，减速组的中位生存期为 6,462 天，而加速组为 4,625 天。

CpGPT 能够强有力地预测个体死亡风险，并有效区分高、低生存概率的人群。

发病率预测：

在训练队列2中，CpGPT 在基线和4年随访时对 9 种疾病（如痴呆、心血管疾病、癌症、糖尿病等）表现出预测能力，其中痴呆和心血管疾病的AUC最高（基线 0.66，随访 0.70）。
在测试队列中，模型在基线和2年随访时对 8 种疾病进行了预测，心绞痛和中风的AUC最高（达 0.74-0.75）。
模型评分还与多种功能性指标显著相关，例如慢性病总数、行动困难、身体脂肪百分比（DEXA扫描）、健康饮食指数等。较高的CpGPT评分与更高的疾病负担和更差的功能状态相关。
模型还与抑郁（CES-D评分）等心理社会指标相关联。

CpGPT 不仅能够预测广泛的衰老相关疾病，还能反映个体的功能状态、代谢健康甚至心理社会健康状况，展示了其作为综合性健康风险评估工具的潜力。