Skip to contents

Colocalization Pipeline (coloc.bf_bf over paired LBF matrices)

Source: R/colocPipeline.R
colocPipeline.Rd

Per-region pipeline that pairs two fine-mapping result collections and runs coloc::coloc.bf_bf per (first-side tuple, second-side tuple) pair to produce per-pair colocalization posterior probabilities PP.H0-PP.H4.

Either side may be a QtlFineMappingResult or a GwasFineMappingResult: QTL-GWAS, QTL-QTL (two molecular phenotypes) and GWAS-GWAS (two diseases) all run through one code path, since nothing below the identity tuple depends on which flavour a side is. The second side may also be handed in as summary statistics (QtlSumStats or GwasSumStats), which are fine-mapped inline.

The argument names keep the QTL / GWAS wording of the common case. What they mean generally is: qtlFineMappingResult is the side whose identity is reported in the unprefixed study / context / trait / method columns, and gwasInput is the side reported in the gwas-prefixed ones. A GWAS side has no context or trait axis, so those two columns are NA for it.

Usage

colocPipeline(
  qtlFineMappingResult,
  gwasInput,
  priors = ColocPriorParam(),
  colocLbfFilterParam = ColocLbfFilterParam(),
  gwasFineMappingParam = GwasFineMappingParam(),
  returnGwasFineMapping = FALSE,
  enrichment = NULL,
  adjustPips = TRUE,
  alleleFlip = TRUE,
  methodArgs = ColocOptions()
)

Arguments

qtlFineMappingResult

The first side: a QtlFineMappingResult or a GwasFineMappingResult (required).

gwasInput

The second side: a QtlFineMappingResult, a GwasFineMappingResult, or the summary statistics to fine-map inline (QtlSumStats or GwasSumStats).

priors

Per-variant prior probabilities, built with ColocPriorParam: p1 (QTL signal), p2 (GWAS signal) and p12 (shared). Mutually exclusive with enrichment, which derives the same three.

colocLbfFilterParam

Credible-set filtering applied before scoring, built with ColocLbfFilterParam. filterLbfCs keeps only effects that produced a credible set; supplying secondary instead runs a concentration filter at that coverage, where an effect is kept only if its credible set spans fewer than nVariants * secondary * concentration variants. Diffuse effects are dropped before the LBF matrix reaches coloc::coloc.bf_bf. priorTol is the default filter's own cutoff (effects with V <= priorTol are dropped) and does not apply once filterLbfCs or secondary selects another.

gwasFineMappingParam

How gwasInput is fine-mapped when it is summary statistics rather than a fine-mapping result, built with GwasFineMappingParam: methods, credibleSetParam, susieRssParam, panelFilterParam, initializeWithSusieInf and fitRetention. Inert when gwasInput is already fine-mapped. credibleSetParam is the one to check: its coverage / maxNumSingleEffects decide the sets whose LBFs are scored, so defaults here against a differently-built QTL side compare two different things.

returnGwasFineMapping

Logical. When TRUE, attach the fine-mapping result computed from gwasInput on the returned object as attribute "gwasFineMapping". Default FALSE.

enrichment

Optional data.frame of per-pair enrichment factors with columns gwasStudy, qtlStudy, qtlContext, enrichment, and optionally gwasContext / gwasTrait (which the join uses when present, and which qtlEnrichmentPipeline emits for a QTL outcome side). Output of qtlEnrichmentPipeline. When non-NULL the pair is scored by enloc rather than coloc: the priors come from that table's colocP1 / colocP2 / colocP12 columns, and the report gains RCP and LCP. Pairs without a matching enrichment row fall back to the baseline priors with a warning. Default NULL (baseline coloc).

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.

methodArgs

Additional arguments forwarded to coloc::coloc.bf_bf, built with ColocOptions. The prior probabilities are set through priors instead, since the enrichment adjustment reads them too.

Value

A ColocResult: one element per tested (first-side credible set, second-side 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, gwasContext, gwasTrait, 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 the pair is scored by enloc instead, and four further columns appear: enrichment and p12Used report the per-pair factor and the shared-signal prior enloc derived, while RCP and LCP are enloc's regional and locus colocalisation probabilities. Those two are the same numbers as PP.H4.abf and PP.H3.abf + PP.H4.abf, reported under enloc's names because that is what readers of enloc output expect to find.

Project it with colocPairs (the flat table this pipeline used to return, also available as as.data.frame), colocVariants, colocCredibleSets or colocGenes.

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: colocLbfFilterParam's filterLbfCs keeps only effects that produced a credible set, its secondary keeps effects at a secondary coverage; otherwise the default filter drops effects whose prior variance is below priorTol.

  • Per-tuple LBF reuse: each second-side tuple's LBF matrix is extracted once and scored against every first-side tuple, so the filtering above is applied once per tuple rather than once per pair.

This pipeline preserves all three.

Second-side input dispatch:

LD reference compatibility check

If ldReference(qtlFineMappingResult) is non-NULL, it must be the same LD reference as the one on gwasInput. The two must be the same kind of reference. Two genotype panels must share their samples and allele orientation; two LdMatrix objects must agree on the LD they share; two GenotypeMixture objects must average the same panels with the same weights. The two need NOT carry the same variants — running summaryStatsQc on the two sides separately normally leaves each reference trimmed to its own surviving variants, and LD is looked up per variant, so a partial overlap only warns (once per session). No shared variant at all, or a mismatch of any of those, is a hard error. When the first side's ldReference is NULL (individual-level fit), the validation is skipped on that side and the second side's panel is what the result carries forward.

Examples

data(qtlFineMappingLbfExample)
data(gwasFineMappingLbfExample)
colocPipeline(qtlFineMappingLbfExample,
  gwasInput = gwasFineMappingLbfExample)
#> ColocResult with 50 colocalized pair(s)
#>   studies     : test_study 
#>   paired with : gwas1 
#>   variants    : 10000 across all pairs
#>   max PP.H4   : 0.001139 
# Either side may be a QTL result. Pairing this collection against itself
# colocalizes its two contexts, and reports the second side's context and
# trait in gwasContext / gwasTrait.
res <- colocPipeline(qtlFineMappingLbfExample,
  gwasInput = qtlFineMappingLbfExample)
unique(colocPairs(res)[, c("context", "gwasContext")])
#> # A tibble: 4 × 2
#>   context  gwasContext
#>   <chr>    <chr>      
#> 1 context1 context1   
#> 2 context1 context2   
#> 3 context2 context1   
#> 4 context2 context2