01使用 coloc 进行共定位分析

整理 GWAS 与 QTL 区域汇总数据,使用 coloc.abf 比较单信号,并在多信号场景中结合 SuSiE 与 LD。

2026-04-06
GWASColocalizationcolocSuSiEQTL
本章目录 · 8

两个性状在同一区域都出现关联信号,不代表它们由同一个 causal variant 驱动。它们可能共享变异,也可能只是各自的 causal variants 处于 linkage disequilibrium(LD)中。

coloc 使用区域内的 summary statistics 比较这些解释。常见场景包括 GWAS 与 eQTL、pQTL 或 mQTL 的比较,也可以比较两个疾病或连续性状。

library(coloc)

五个假设

经典 coloc.abf() 为一个区域计算五种假设的 posterior probability:

假设 含义
H0 两个性状在区域内都没有关联
H1 只有性状 1 有关联
H2 只有性状 2 有关联
H3 两个性状都有信号,但由不同 causal variants 驱动
H4 两个性状的信号由同一个 causal variant 驱动

分析重点通常是比较 PP.H3PP.H4,而不是只问 PP.H4 是否超过某个固定阈值。若区域对其中一个性状缺乏足够关联证据,H0、H1 或 H2 可能占主导,此时不能把较低的 H4 简单解释为“不共定位”。

输入数据结构

每个性状分别整理为一个 named list。优先使用效应值与方差:

字段 含义 说明
beta SNP effect estimate varbeta 配套
varbeta effect estimate variance 通常为 se^2
snp variant identifier 两个数据集应匹配同一批变异
position genomic position 检查区域和绘图时使用
type "quant""cc" 定量或病例对照性状
sdY outcome standard deviation 定量性状需要,部分情况下可估计
s case fraction 病例对照性状需要
N sample size 缺少部分参数或使用 p 值输入时需要
MAF minor allele frequency 参数估计及 p 值输入时需要

beta 和标准误时:

dataset1 <- list(
  beta = gwas$beta,
  varbeta = gwas$se^2,
  snp = gwas$snp,
  position = gwas$position,
  type = "cc",
  s = gwas$n_case / gwas$n_total,
  N = gwas$n_total,
  MAF = gwas$maf
)

dataset2 <- list(
  beta = eqtl$beta,
  varbeta = eqtl$se^2,
  snp = eqtl$snp,
  position = eqtl$position,
  type = "quant",
  N = eqtl$n,
  MAF = eqtl$maf,
  sdY = eqtl$sd_y
)

缺少 beta 与 variance 时,coloc 也可以使用 pvaluesNMAF 等信息进行近似,但完整效应估计通常更合适。不要把标准误直接传给 varbeta

varbeta <- se^2

先对齐区域和变异

两个数据集必须使用相同 genome build,并在同一分析区域中取交集:

common_snps <- intersect(gwas$snp, eqtl$snp)

gwas_region <- gwas[gwas$snp %in% common_snps, ]
eqtl_region <- eqtl[eqtl$snp %in% common_snps, ]

随后按同一 SNP 顺序排列,并检查:

  • rsID 或 chr:position:alleles 是否指向同一变异;
  • effect allele 与 other allele 是否已经 harmonize;
  • 是否混用了 genome builds;
  • palindromic variants 的方向能否可靠判断;
  • MAF、样本量和 case fraction 是否来自相应数据集;
  • 区域是否完整覆盖两边的关联信号,而不是只保留 lead SNP。

coloc 使用区域内完整的关联模式。只拿各自显著 SNP 的交集,会破坏这个模式并产生误导。

检查数据集

内置测试数据可以展示最小结构:

data(coloc_test_data)

D1 <- coloc_test_data$D1
D2 <- coloc_test_data$D2

minimum_data <- D1[c(
  "beta", "varbeta", "snp", "position", "type", "sdY"
)]

check_dataset(minimum_data)
plot_dataset(minimum_data)

对真实数据分别运行检查:

check_dataset(dataset1)
check_dataset(dataset2)

检查通过只说明结构满足函数要求,不证明 allele alignment、研究人群或生物学区域选择正确。

单信号分析:coloc.abf()

经典分析假设每个性状在该区域最多有一个 causal signal:

result <- coloc.abf(
  dataset1 = dataset1,
  dataset2 = dataset2
)

result$summary

查看对 H4 贡献较高的 variants:

head(
  result$results[order(result$results$SNP.PP.H4, decreasing = TRUE), ],
  10
)

SNP.PP.H4 是在 H4 成立条件下各 SNP 成为共享 causal variant 的相对支持,不应脱离区域级 PP.H4 单独报告。

Prior sensitivity

H4 的 posterior probability 会受 priors 影响,尤其是 p12,即一个 variant 同时关联两个性状的先验概率。

result <- coloc.abf(
  dataset1 = dataset1,
  dataset2 = dataset2,
  p12 = 1e-6
)

sensitivity(result, rule = "H4 > 0.5")

规则应在分析计划中说明,并结合 H3、区域信号强度和研究目的解释。若结论只在很窄的 prior 范围内成立,应将这种不稳定性报告出来。

多信号区域:coloc.susie()

同一区域可能存在多个独立信号。此时单 causal variant 假设不合适,可以提供 LD matrix,先用 SuSiE 分解信号,再进行共定位。

check_dataset(dataset1, req = "LD")
check_dataset(dataset2, req = "LD")

susie1 <- runsusie(dataset1)
susie2 <- runsusie(dataset2)

result_susie <- coloc.susie(susie1, susie2)
result_susie$summary

LD matrix 必须:

  • 为 SNP × SNP 方阵;
  • 行名、列名和 dataset 中的 snp 完全对应;
  • SNP 顺序一致;
  • 来自与 summary statistics 祖源尽量匹配的参考人群。
identical(rownames(dataset1$LD), dataset1$snp)
identical(colnames(dataset1$LD), dataset1$snp)

参考 LD 不匹配会改变信号分解。coloc.susie() 能处理多个信号,但不能修复错误的 LD 或 allele alignment。

结果边界

共定位支持两个关联信号与共享 causal variant 相容,但它本身不证明:

  • 一个性状导致另一个性状;
  • 对应基因就是疾病机制中的 effector gene;
  • 共享变异只通过所测分子性状发挥作用;
  • 不同组织、细胞类型或祖源中的信号相同。

报告时应写清楚区域、build、样本量、性状类型、variant 数量、priors、LD 来源、H0–H4 posterior probabilities,以及是否采用多信号模型。

参考:coloc documentationsusieR documentation