Causal TWAS Pipeline (cTWAS, multi LD block)

Pipeline that hands a per-block set of GwasSumStats of GWAS Z-scores together with the matching per-block per-gene TWAS weights and LD sketches to ctwas::ctwas_sumstats, producing per-gene posterior inclusion probabilities for causal genes. Optionally accepts a precomputed TWAS-Z GRanges from causalInferencePipeline as the z_gene input so the per-gene Z is not recomputed inside ctwas.

Usage

ctwasPipeline(
  gwasSumStats,
  twasWeights,
  twasZ = NULL,
  fineMappingResult = NULL,
  method = NULL,
  thin = 0.1,
  niterPrefit = 3L,
  niter = 30L,
  L = 5L,
  groupPriorVarStructure = c("shared_type", "shared_context", "shared_nonSNP",
    "shared_all", "independent"),
  ncore = 1L,
  twasWeightCutoff = 0,
  csMinCor = 0.8,
  minPipCutoff = 0,
  maxNumVariants = Inf,
  fallbackToPrefit = FALSE,
  keepSnps = FALSE,
  mergeBoundary = FALSE,
  mergePipThresh = 0.5,
  mergeFilterCs = FALSE,
  mergeMaxSNP = Inf,
  ...
)

Arguments

gwasSumStats: NAMED LIST of GwasSumStats keyed by region_id (at least two entries). Each must have getQcInfo() non-empty.
twasWeights: The per-gene weight source. Either (a) a FLAT TwasWeights / QtlFineMappingResult (or a homogeneous list of them) carrying region provenance — each gene is placed into its home LD block internally by start(region) (matching cTWAS's p0 assignment rule); or (b) a pre-bucketed NAMED LIST keyed by region_id (keys a SUBSET of gwasSumStats's), used as-is. Blocks without any TWAS weights still contribute their SNP-level signal to ctwas's joint group prior estimate (the legacy whole-chromosome pattern where only a few of many LD blocks carry gene weights). A gene whose cis span straddles a block boundary is homed by its single anchor; the cross-block signal is cTWAS's boundary-gene concern (mergeCtwasBoundaryRegions), not placement.
twasZ: Optional GRanges of TWAS Z-scores (output of causalInferencePipeline). When supplied, the per-(trait, context) Z is used as the z_gene input to ctwas_sumstats so it is not recomputed.
fineMappingResult: Optional QtlFineMappingResult or GwasFineMappingResult carrying the per-variant PIP and credible-set membership data used by the CS / PIP rescue filters (csMinCor and minPipCutoff). When NULL (default) the smart filters are no-ops; only the magnitude filter (twasWeightCutoff) and the per-gene cap (maxNumVariants, ordered by |weight|) apply.
method: Optional character (length 1). Picks which TWAS method's weights to feed into ctwas for each (study, context, trait) gene. When NULL (default): use "ensemble" if that method is present across the weight sources; otherwise use the sole method when only one is present; otherwise run every method as an independent cTWAS run (one CtwasResult row-set per method). Passing the name explicitly (e.g. "mrash") restricts the run to that single method.
thin, niterPrefit, niter, L: Pass-throughs to ctwas::ctwas_sumstats.
groupPriorVarStructure: Pass-through (defaults "shared_type").
ncore: Number of cores. Default 1.
twasWeightCutoff: Numeric (length 1). Drop variants with |weight| < twasWeightCutoff from each gene's weight matrix before ctwas sees it. Default 0 (no filter).
csMinCor: Numeric (length 1). When fineMappingResult is provided, variants belonging to any 95% credible set with purity (min_abs_corr) >= csMinCor are marked as must-keep and survive the per-gene cap. Default 0.8. Ignored without a fineMappingResult.
minPipCutoff: Numeric (length 1). When fineMappingResult is provided, variants with PIP greater than minPipCutoff are marked as must-keep and survive the per-gene cap. Default 0 (no PIP rescue). Ignored without a fineMappingResult.
maxNumVariants: Numeric (length 1). Cap on per-gene variant count. When the gene has more variants than this, keep all must-keep variants and fill remaining slots by descending PIP (when available) or descending |weight|. Default Inf (no cap).
fallbackToPrefit: Logical (length 1). Forwarded to estCtwasParam. When TRUE, ctwas's accurate-EM NaN failure is recovered by falling back to the prefit estimates (mirrors the legacy ctwas_2 workaround on underpowered data). Default FALSE.
keepSnps: Logical (length 1). When TRUE, retain the context-agnostic SNP background of each run as one extra CtwasResult row (study = context = "SNP", mirroring cTWAS's own "SNP" group) so the full ctwas output is reconstructable from getFinemap / getSusieAlpha. Default FALSE — the SNP rows are the null background and are dropped from the structured gene-level result.
mergeBoundary: Logical (length 1). When TRUE, run mergeCtwasBoundaryRegions after fine-mapping each run: a high-PIP gene whose cis window straddles an LD-block boundary has its adjacent regions merged and re-fine-mapped (the legacy default-off ctwas_3 post-processing). Default FALSE.
mergePipThresh: Numeric (length 1). PIP threshold for selecting which boundary genes to merge (mergeCtwasBoundaryRegions pipThresh). Default 0.5. Ignored unless mergeBoundary = TRUE.
mergeFilterCs: Logical (length 1). Require the boundary gene to be in a credible set to be selected. Default FALSE. Ignored unless mergeBoundary = TRUE.
mergeMaxSNP: Numeric (length 1). Per-merged-region SNP cap. Default Inf. Ignored unless mergeBoundary = TRUE.
...: Additional arguments forwarded to ctwas::ctwas_sumstats.

Value

A CtwasResult collection: one row per (gwasStudy, study, context, method). A single-context run is one row per method; a multi-context (joint) run emits per-context rows sharing the same jointContexts set and the jointly-estimated group priors. Each row's CtwasResultEntry payload carries that context's per-gene fine-mapping posteriors (finemap), the run's param, and its regionInfo. For the raw ctwas::finemap_regions list (e.g. to feed mergeCtwasBoundaryRegions), call the granular assembleCtwasInputs \(\to\) estCtwasParam \(\to\) screenCtwasRegions \(\to\) finemapCtwasRegions path instead.

LD block convention

Inputs are NAMED LISTS keyed by region_id (list(block1 = gss1, block2 = gss2, ...)). Per-block region_info, LD_map, and snp_map entries are built automatically from each block's LD sketch and concatenated before the call to ctwas::ctwas_sumstats. A single-block input is rejected: cTWAS's EM cannot converge on a single region, so callers must supply at least two blocks.

LD-sketch identity check

Per block: getLdSketch(twasWeights) (when non-NULL) must match getLdSketch(gwasSumStats). Mismatch is a hard error.