Colocalization Pipeline (coloc.bf_bf over QTL + GWAS LBF matrices)

Per-region pipeline that pairs a QTL QtlFineMappingResult with a GWAS fine-mapping result (either supplied directly as a GwasFineMappingResult or computed inline from a GwasSumStats) and runs coloc::coloc.bf_bf per (QTL tuple, GWAS tuple) pair to produce per-pair colocalization posterior probabilities PP.H0-PP.H4.

Usage

colocPipeline(
  qtlFineMappingResult,
  gwasInput,
  filterLbfCs = FALSE,
  filterLbfCsSecondary = NULL,
  filterLbfCsConcentration = 0.5,
  priorTol = 1e-09,
  p1 = 1e-04,
  p2 = 1e-04,
  p12 = 5e-06,
  finemappingMethods = "susie",
  returnGwasFineMapping = FALSE,
  enrichment = NULL,
  p12Max = 0.001,
  adjustPips = TRUE,
  alleleFlip = TRUE,
  ...
)

Arguments

qtlFineMappingResult

A QtlFineMappingResult (required).

gwasInput

Either a GwasSumStats or a GwasFineMappingResult.

filterLbfCs

Logical. When TRUE (and filterLbfCsSecondary is NULL), keep only effects that produced a credible set (trimmedFit$sets$cs_index). Default FALSE.

filterLbfCsSecondary

Optional secondary coverage (numeric in $(0, 1)$). When supplied, run a credible-set concentration filter at this coverage level instead of filterLbfCs: each L-effect's credible set must span fewer than nVariants * filterLbfCsSecondary * filterLbfCsConcentration variants to be kept. Effects with diffuse credible sets are dropped before the LBF matrix is passed to coloc::coloc.bf_bf. Overrides filterLbfCs when set.

filterLbfCsConcentration

Numeric in $(0, 1)$; the concentration factor in the cutoff above. With the default 0.5 a 50% credible set is kept only if it spans fewer than 25% of the locus's variants. Only consulted when filterLbfCsSecondary is non-NULL. Default 0.5.

priorTol

Prior-variance cutoff for the default filter: effects with V <= priorTol are dropped. Ignored when either filterLbfCs or filterLbfCsSecondary is in use. Default 1e-9.

p1

Prior probability of QTL signal per variant. Default 1e-4.

p2

Prior probability of GWAS signal per variant. Default 1e-4.

p12

Prior probability of shared signal per variant. Default 5e-6.

finemappingMethods

Character vector forwarded to fineMappingPipeline when gwasInput is a GwasSumStats. Default "susie".

returnGwasFineMapping

Logical. When TRUE, attach the computed GwasFineMappingResult on the returned data frame as attribute "gwasFineMapping". Default FALSE.

enrichment

Optional data.frame of per-(gwasStudy, qtlStudy, qtlContext) enrichment factors with columns gwasStudy, qtlStudy, qtlContext, enrichment. Output of qtlEnrichmentPipeline. When non-NULL, each pair's p12 prior is scaled to min(p12 * (1 + enrichment), p12Max) (the enrichment-informed colocalization variant, "enloc"). Pairs without a matching enrichment row fall back to the baseline p12 with a warning. Default NULL (baseline coloc).

p12Max

Numeric scalar. Maximum value for the enrichment-adjusted p12 prior. Default 1e-3. Ignored when enrichment = NULL.

adjustPips

Logical, default TRUE. When TRUE, before any per-pair inference the QTL and GWAS fine-mapping result collections are reconciled with intersectVariants, so both sides are scored on the variants they share and their PIPs are renormalized to that shared set. This matters in two scenarios: (1) the user declined to impute missing variants in the GWAS SumStats and the QTL fine-mapping input has additional variants; (2) the GWAS fine-mapping result contains variants not present in the QTL fine-mapping result. Pass FALSE to use the FMRs as supplied.

How much of each effect's posterior survived that restriction is reported per result row as qtlRetainedMass / gwasRetainedMass; a low value means the effect was largely built on variants the other side does not carry, so its colocalization evidence rests on little retained signal. Both are NA when adjustPips = FALSE, since nothing was measured.

alleleFlip

Logical, default TRUE. When TRUE, align LBF columns between the QTL and GWAS by (chrom, pos) with ref/alt swaps recognized (LBF is coding-invariant, so no sign change is needed); when FALSE, match on exact alleles only, so a ref/alt swap is treated as a distinct variant.

...

Additional arguments forwarded to coloc::coloc.bf_bf.

Value

A ColocResult: one element per tested (QTL credible set, GWAS credible set, block) pair, holding that pair's aligned variants with their SNP.PP.H4. Pair-level metadata carries the identity columns (study, context, trait, method, gwasStudy, gwasMethod), the block and stable credible-set ids (blockId, qtlCs, gwasCs), the standard coloc fields (idx1, idx2, nSnps, hit1, hit2, PP.H0.abf ... PP.H4.abf) and the reconciliation diagnostics qtlRetainedMass / gwasRetainedMass. When enrichment is supplied, two additional columns enrichment and p12Used report the per-pair factor and the prior actually passed to coloc::coloc.bf_bf.

Project it with getColocPairs (the flat table this pipeline used to return, also available as as.data.frame), getColocVariants, getColocCredibleSets or getColocGenes.

Why `coloc.bf_bf` and not `coloc.susie`

The prior colocWrapper (now stubbed) used coloc::coloc.bf_bf on the SuSiE lbf_variable matrices directly. That choice carries three behaviours that coloc::coloc.susie does not expose:

fSuSiE support: the LBF matrix lives at a different slot for fSuSiE fits (fsusie_result$lBF) and gets stacked into a single matrix.
Effect filtering: filterLbfCs keeps only effects that produced a credible set; filterLbfCsSecondary keeps effects at a secondary coverage; otherwise the default filter drops effects whose prior variance is below priorTol.
Multiple-GWAS batching: when several GWAS fine-mapping rows fall in the same region they are merged into one combined LBF matrix per QTL pair (one coloc.bf_bf call covers them all).

This pipeline preserves all three.

GWAS input dispatch:

gwasInput is a GwasSumStats: GWAS fine-mapping is performed inline by fineMappingPipeline with the supplied finemappingMethods (default "susie").
gwasInput is a GwasFineMappingResult: used directly; no inline fine-mapping.

LD-sketch identity check

If getLdSketch(qtlFineMappingResult) is non-NULL, it must match the LD sketch on gwasInput. Mismatch is a hard error. When the QTL FMR's ldSketch is NULL (individual-level fit), the validation is skipped on the QTL side.

Examples

data(qtlFineMappingLbfExample)
data(gwasFineMappingLbfExample)
colocPipeline(qtlFineMappingLbfExample,
  gwasInput = gwasFineMappingLbfExample)
#> ColocResult with 50 colocalized pair(s)
#>   QTL studies : test_study 
#>   GWAS studies: gwas1 
#>   variants    : 10000 across all pairs
#>   max PP.H4   : 0.001139