08 表观组与转录组的整合
表观修饰（开放染色质、TF 结合、甲基化）最终要通过影响基因表达来发挥功能。把表观组数据和转录组数据放在一起看，能回答"哪些表观变化真正驱动了表达变化"。
常见整合策略
| 策略 | 输入 | 输出 | 工具 |
|---|---|---|---|
| Peak-gene 关联 | 差异 peak + 差异基因 | 重叠基因列表 | ChIPseeker + 自定义脚本 |
| 相关性分析 | peak 信号矩阵 + 表达矩阵 | peak-gene 相关性 | cor.test / LOLA |
| 调控网络推断 | motif + 表达 + peak | TF → target 网络 | SCENIC / pySCENIC |
| 多组学因子分析 | 多层矩阵 | 共变因子 | MOFA2 / mixOmics |
Peak-gene 关联（最简单）
library(ChIPseeker)
# 差异 peak 注释到最近基因
db_anno <- annotatePeak(db_peaks, TxDb = txdb)
db_genes <- unique(as.data.frame(db_anno)$SYMBOL)
# 差异基因（来自 DESeq2）
de_genes <- res_df$SYMBOL[res_df$padj < 0.05]
# 交集
overlap <- intersect(db_genes, de_genes)
cat("Overlap:", length(overlap), "genes\n")
# Fisher 检验看是否显著富集
fisher.test(matrix(c(
length(overlap),
length(setdiff(db_genes, de_genes)),
length(setdiff(de_genes, db_genes)),
total_genes - length(union(db_genes, de_genes))
), nrow = 2))
如果 overlap 显著大于随机期望，说明表观变化和表达变化确实有关联。
方向一致性检查
更严格的验证：不仅看"有没有重叠"，还看"方向是否一致"：
# 合并 peak FC 和 gene FC
merged <- inner_join(
data.frame(gene = db_genes_df$SYMBOL, peak_fc = db_genes_df$Fold),
data.frame(gene = de_df$SYMBOL, rna_fc = de_df$log2FoldChange)
)
# 散点图：peak FC vs RNA FC
ggplot(merged, aes(x = peak_fc, y = rna_fc)) +
geom_point() +
geom_smooth(method = "lm") +
labs(x = "Peak log2FC (ChIP/ATAC)", y = "RNA log2FC")
正相关 = 开放区域增加的基因表达也增加（符合预期）。如果是甲基化数据，启动子区域应该是负相关（甲基化增加 → 表达下降）。
SCENIC：从 scATAC + scRNA 推断调控网络
如果有配对的单细胞数据（10x Multiome 或分别测的 scRNA + scATAC），SCENIC+ 能推断出"哪个 TF 通过哪个增强子调控哪个基因"：
# Python (pySCENIC+)
import scenicplus
# 输入：scRNA AnnData + scATAC AnnData + motif 数据库
# 输出：TF → enhancer → gene 的三元组网络
这是目前单细胞表观组整合的最前沿方向，计算量大但信息量也最大。
实用建议
- 先做简单的 peak-gene overlap，确认方向一致性
- 如果 overlap 显著且方向一致，再做更复杂的网络推断
- 多组学整合的结果要用独立实验验证（比如 CRISPRi 敲掉某个增强子看表达是否下降）
- 不要过度解读"相关性 = 因果性"
整合分析的预期值要校准
新手做 ATAC + RNA 整合常常失望："只有 30% 的差异 ATAC peak 关联的基因在 RNA-seq 里也是差异的"。这其实是正常的，甚至是好的：
| 重叠比例 | 解读 |
|---|---|
| > 50% | 实验设计强相关（同一组织、同一时间点） |
| 20-40% | 正常，差异 ATAC ≠ 一定改变表达 |
| 5-15% | 偏低，但生物学合理（很多 ATAC 变化是 priming，表达变化滞后） |
| < 5% | 真的有问题（条件混淆、批次没去干净、差异分析阈值不合理） |
ATAC 变化先于 RNA 变化、增强子调控不一定立即翻译成表达 — 这些都是真实的生物学。期望"100% 一致"是统计学上的过度自信。
常见坑
坑 1：peak-gene 关联用 nearest-gene 当唯一标准
最近基因 ≠ 真调控基因。严格的关联需要 Hi-C / ABC model / CRISPRi 验证。教程级用 nearest 可以，发表级要补 chromatin loop 数据或至少加距离过滤（比如 < 100kb）。
坑 2：不区分启动子 vs 增强子 peak
启动子 peak 和增强子 peak 的 RNA 关联机制不同（启动子直接，增强子要先经过 looping）。整合时应该按 peak 类型分开做关联，混着看会稀释信号。
坑 3：方向一致性看的是 raw correlation
raw 表达矩阵和 raw peak signal 的 correlation 受 batch / library size 影响大。要先 normalize（vst / TMM）再算 correlation。
坑 4：把"显著富集"误读成"重要驱动"
Fisher test 显著只说明重叠不是随机的，不说明这个 overlap 集合里的基因就是关键调控者。后续要做 motif 富集 + KO 验证才能锁定 driver。
坑 5：SCENIC / MOFA2 跑出结果不验证就上文章
这类高级方法对参数敏感，跑两次结果会不一样。任何整合分析的结论都要用独立数据集（公开数据 + own 验证）confirm，不要依赖单次跑出的网络。
下一步
接着深入：
- 多组学整合 overview — 系统的多组学方法学专栏
- 多组学整合 module01 — 项目设计与样本匹配
横向延伸：
- 03 DiffBind 差异结合分析 — 整合的输入是各 omics 的差异表
- bulk RNA-seq 02 DESeq2 — 表达层的差异分析
参考资源
离线资料下载
手册 HTML / PDF 已在后台预生成，点击后直接下载网站静态资源。