Bio Causal Genomics Colocalization Analysis

Version Compatibility

Reference examples tested with: coloc 5.2.3+, susieR 0.12.35+, hyprcoloc 1.0+ (GitHub jrs95/hyprcoloc), SMR 1.3.1+ (CLI, cnsgenomics.com), eCAVIAR 2.2+ (compiled from caviar/eCAVIAR repo), PWCoCo 1.0+ (jwr-git/pwcoco), moloc 0.1+ (clagiamba/moloc), SharePro_coloc 7.0+ (zhwm/SharePro_coloc), R >= 4.1.

Before using code patterns, verify installed versions match. If versions differ:

• R: packageVersion('coloc'); check ?coloc.abf, ?coloc.susie, ?runsusie
• CLI: smr --version, pwcoco --help, sharepro_coloc.py --help

If code throws AttributeError, NULL list elements, or Error in coloc.abf: dataset must have..., introspect the installed package signature and adapt the example rather than retrying.

Colocalization Analysis

"Test whether my GWAS signal and an eQTL share the same causal variant" -> Compute Bayesian posterior probabilities over five hypotheses (H0 neither, H1 trait-1-only, H2 trait-2-only, H3 distinct causal variants, H4 shared causal variant) to discriminate true causal overlap from LD-driven coincidence, then run sensitivity analysis over the p12 prior.

• R (single-causal, fastest): coloc::coloc.abf(dataset1, dataset2, p12=5e-6) -> coloc::sensitivity(res, 'H4 > 0.75')
• R (multi-causal, needs LD): runsusie(d1) -> runsusie(d2) -> coloc.susie(s1, s2) -> per-credible-set PP
• R (many traits, single-causal cluster): hyprcoloc::hyprcoloc(effect.est = betas_mat, effect.se = ses_mat, trait.names = ..., snp.id = ...) -> trait clusters
• CLI (causality vs linkage): smr --bfile ref --gwas-summary g.ma --beqtl-summary eqtl.besd --out smr -> SMR p + HEIDI p
• CLI (allelic heterogeneity): eCAVIAR eCAVIAR -l ld1 -l ld2 -z z1 -z z2 -o out -c 2 -> CLPP per SNP
• CLI (conditional): PWCoCo conditions on each independent signal via GCTA-COJO then runs pairwise coloc.abf

Algorithmic Taxonomy

Methodology evolves; verify the current Open Targets Genetics, eQTL Catalogue, and FinnGen colocalization pipelines before locking parameters. Open Targets uses coloc.abf at PP.H4 >= 0.75 with p12 = 1e-5; FinnGen uses coloc.susie at PP.H4 >= 0.8 with in-sample LD.

Decision Tree by Scenario

Per-Method Failure Modes

coloc.abf -- PP.H3 inflation under multiple causal variants

Trigger: Locus has 2+ independent causal signals in moderate LD (r2 ~ 0.3-0.6).

Mechanism: The single-causal-variant assumption forces the model to allocate posterior mass to H3 (distinct causal variants) whenever the per-SNP Bayes factors for the two top SNPs do not align.

Symptom: Visual co-localization in LocusZoom looks convincing, but result$summary['PP.H3.abf'] dominates over PP.H4; sensitivity() shows PP.H4 stays low across the entire p12 grid.

Fix: Run coloc.susie (or eCAVIAR or PWCoCo) to allow multiple causal variants. If coloc.susie returns multiple credible sets with one pair showing PP.H4 > 0.75, this is real allelic heterogeneity not failure.

coloc.susie -- LD reference mismatch

Trigger: Z-scores from GWAS / eQTL of one ancestry, LD matrix from 1000 Genomes EUR (or any non-matched reference).

Mechanism: SuSiE assumes z-scores and the supplied LD are jointly consistent. Ancestry mismatch or sample-size mismatch produces a non-positive-definite implicit covariance; SuSiE responds by returning spurious credible sets that include LD-mismatched SNPs.

Symptom: susieR::estimate_s_rss(z, R, n) returns lambda > 0.05; susieR::kriging_rss flags off-diagonal SNPs with extreme studentized residuals; credible sets are oddly large (50+ SNPs) or include SNPs distant in LD from the lead.

Fix: Use in-sample LD when at all possible (per-cohort plink --r square). If reference must be external, match ancestry (1KG superpopulation) and superpopulation-stratify. Run estimate_s_rss and report lambda; if > 0.05, drop the locus or switch to coloc.abf.

coloc default p12 too liberal for trans-eQTL

Trigger: Applying p12 = 1e-5 (the default) to a trans-eQTL / GWAS pair.

Mechanism: The default p12 was calibrated for cis-eQTL where biological proximity makes shared causality reasonable. For trans-eQTL, prior probability of shared causality is much lower; uniform p12 over-favours H4.

Symptom: PP.H4 > 0.8 reported, but sensitivity() reveals PP.H4 falls below 0.5 for p12 < 1e-5; replication in independent data fails.

Fix: Operational definition: "trans" = >5 Mb from TSS or different chromosome. Default p12=1e-5 over-favors H4 for trans (genome-rare biology). For trans: lower p12 to 5e-6 or 1e-6 AND raise PP.H4 threshold to >= 0.8 (compensate for higher FP risk). Cross-reference Vosa 2021 Nat Genet 53:1300 (eQTLGen trans) for empirical patterns.

MHC / HLA + chr 8 inversion -- single-causal assumption breaks

Trigger: Locus within chr6:25-35 Mb (extended MHC, hg38), or chr8:8.1-11.9 Mb (chr 8 inversion, hg38).

Mechanism: The MHC contains classical HLA genes with extreme long-range LD (r2 > 0.5 over many Mb), multiple independent causal haplotypes, and structural variation. The chr 8p23.1 inversion similarly produces long-range LD across megabases of polymorphic inversion alleles. The single-causal-variant assumption is biologically wrong in both regions.

Symptom: coloc.abf almost always returns PP.H3 or fragmented PP across H1/H2/H3/H4 even when the underlying biology is well-established (e.g. HLA-DRB1 in autoimmune GWAS).

Fix (MHC): Use HLA-imputed classical alleles via SNP2HLA / HIBAG / HLA-TAPAS, then HLA-coloc (Butler-Laporte 2024 medRxiv) -- NOT coloc on SNPs in MHC. OR exclude MHC from genome-wide coloc and report HLA association at the haplotype/allele level. Fix (chr 8 inversion): Exclude chr8:8.1-11.9 Mb or pre-condition on inversion genotype before coloc. Never report a single coloc PP.H4 in either region without this caveat.

Lead-SNP swap and window bias

Trigger: The two traits have different lead SNPs at the same locus; analyst centres each window on the trait-specific lead.

Mechanism: coloc PP is sensitive to the SNPs in the window; centring on different leads gives different per-SNP overlap and biases toward H3.

Symptom: Re-centring the window on the GWAS lead vs the eQTL lead produces qualitatively different PP.H4.

Fix: Use a SINGLE window (typically +/- 500 kb or 1 Mb) centred on the joint top-variant (the SNP with the lowest min-p across both traits), or on the GWAS lead consistently. Report PP under multiple centring choices; flag the locus if PP swings > 0.2 across centrings.

Underpowered eQTL (N < 200)

Trigger: Small eQTL discovery (e.g. tissue-specific bulk study, N < 200; per-cell-type sc-eQTL).

Mechanism: With low N, varbeta is large; the eQTL's per-SNP Bayes factors are flat; the joint likelihood concentrates on H0 or H1 (GWAS-only).

Symptom: PP.H0 or PP.H1 dominates; the eQTL panel shows visible signal but coloc cannot resolve causal vs noise.

Fix: Use eQTLGen (N ~ 31k whole-blood) or GTEx v8 (N ~ 70-700 per tissue) where possible. For rare cell types, accept the limitation and report the locus as underpowered rather than claim absence of colocalization.

Reference QTL panel choice

GTEx v8 (838 donors, 49 tissues, 2020) is the current PredictDB-supported standard. GTEx v10 (released 2024) has limited harmonisation and is not yet PredictDB-default. eQTLGen blood meta-eQTL (N ~ 31k) wins on sample size for blood cis-eQTL discovery, beating any single tissue on power. Always pin version in methods (e.g. "GTEx v8 MASHR-EUR, PredictDB release 2022-01").

PP.H4 Threshold Framework

Operational rule: Three operational tiers map onto the most common downstream uses: (a) >= 0.7 when PP.H4 is one of several lines of triangulating evidence (TWAS + coloc, cis-MR + coloc, effector-gene multi-evidence) -- this is the threshold downstream skills (causal-genomics/transcriptome-wide-association, causal-genomics/mendelian-randomization cis-MR, causal-genomics/effector-gene-prioritization, causal-genomics/proteome-mr-drug-target) require; (b) >= 0.8 for standard peer-reviewed publication as a stand-alone coloc claim (Wallace 2020); (c) >= 0.95 for industry / clinical drug-target submission. Open Targets and FinnGen pipelines screen at >= 0.75 but downstream publication-grade coloc claims should clear >= 0.8 and triangulation claims >= 0.7. ALWAYS report PP.H3 alongside PP.H4 -- a locus with PP.H4 = 0.6, PP.H3 = 0.3 is qualitatively different from PP.H4 = 0.6, PP.H3 = 0.05 (the former is real ambiguity over single vs distinct causal; the latter is underpowered evidence). Run coloc::sensitivity() and report the p12 range over which PP.H4 stays above the threshold.

Default Priors and the p12 Sensitivity Question

The p12/p1 ratio (= 0.1 under defaults) is the prior odds of colocalization given a trait-1 association. Wallace 2020 (PLoS Genet 16:e1008720) showed default p12 = 1e-5 is too liberal for many real-world settings and recommended sensitivity analysis as standard practice. Pullin & Wallace 2025 (PLoS Genet 21:e1011697) extended this with variant-specific priors weighted by functional annotation.

p12 Sensitivity Grid

Operational rule: Require PP.H4 to remain above threshold across at least 3 adjacent grid points; report the lowest p12 at which PP.H4 >= 0.75. Use coloc::sensitivity(result, rule = 'H4 > 0.75') for the diagnostic plot.

Required reporting: PP.H4 at default priors + p12 range over which PP.H4 stays above threshold.

eCAVIAR CLPP Threshold Framework

CLPP (Colocalization Posterior Probability) is the per-SNP product of the two per-trait fine-mapping posteriors. Threshold conventions:

• Hormozdiari 2016 AJHG 99:1245 used CLPP >= 0.01 (validated against null simulations).
• 2024 GTEx / Open Targets pipelines use CLPP >= 0.05.
• High-confidence claims require CLPP >= 0.1.
• Report both sum-CLPP across the credible set AND max-CLPP at any single SNP -- the two answer different questions (locus-level vs lead-SNP-level confidence).

eCAVIAR -l ld_gwas.ld -l ld_eqtl.ld \
        -z gwas.z -z eqtl.z \
        -o coloc_out -c 2     # -c = max independent causal variants per trait
# Output: per-SNP CLPP in coloc_out_col file; report sum and max

LD Matrix Construction for coloc.susie

Requirements:

• Signed Pearson r (not r2). coloc.susie expects directional LD; squared LD silently inverts effect-direction inference.
• Ancestry-matched to GWAS / eQTL ancestry. EUR LD on AFR z-scores produces spurious credible sets.
• SNP-order-aligned to the beta vector and named to match (row/column names = SNP IDs).
• Positive semi-definite. Numerical-noise negative eigenvalues must be repaired.
• Effective N sample-size-matched to the trait being fine-mapped (provide via runsusie(..., n = N)).

# plink2 phased r (signed Pearson); square matrix output
plink2 --pfile 1KG_EUR \
    --extract snps.txt --chr 6 --from-bp X --to-bp Y \
    --r-phased square --out locus_ld

# Alternative: in-sample LD from BED via bigsnpr
R <- bigsnpr::snp_cor(snp_obj$genotypes, ind.col = locus_snps)
# PSD repair if negative eigenvalues from numerical noise
R <- as.matrix(Matrix::nearPD(R)$mat)
dimnames(R) <- list(snp_ids, snp_ids)

Critical: Row and column order of R MUST match SNP order in the beta vector -- silent failure otherwise. The SuSiE objective stays finite under mis-ordering and returns nonsense credible sets. Verify with stopifnot(rownames(R) == names(beta)) before runsusie. Cross-reference causal-genomics/fine-mapping for the full LD diagnostic block (estimate_s_rss lambda < 0.05, kriging_rss outlier inspection).

SMR vs coloc Reconciliation

SMR (Zhu 2016) and coloc test related but non-identical questions:

• SMR tests pleiotropy vs linkage: does the top eQTL SNP show a GWAS effect explainable by its eQTL effect (pleiotropic / causal) or does the GWAS effect come from a different SNP in LD (linkage)?
• coloc tests shared vs distinct causal variants over an entire window of SNPs.
• HEIDI is SMR's heterogeneity test; null hypothesis is single shared causal SNP. Zhu 2016 Nat Genet 48:481 specifies HEIDI p > 0.05 (NOT 0.01) as non-rejection of single shared causal. HEIDI p > 0.05 does NOT prove shared causality -- only that data cannot reject it; pair with SMR p Bonferroni-corrected across probes. When LD between causal SNPs > 0.7, HEIDI loses power same as coloc.

When LD between two true causal SNPs is high (r2 > 0.7), both SMR/HEIDI and coloc.abf lose discriminatory power: SMR cannot pick which of the LD-tied SNPs is causal, and coloc.abf cannot reject H4 even if biology is two-distinct-causal. coloc.susie + ancestry-matched LD is the modern resolution.

Operational rule: SMR + HEIDI is appropriate when the question is "does this eQTL gene mediate the GWAS effect at all?" coloc is appropriate when the question is "do the two traits share a causal variant in this window?" Run both; agreement (significant SMR + non-rejected HEIDI + PP.H4 >= 0.75) is high-confidence; disagreement requires inspection (often the multi-causal / LD scenario above).