Causal Inference Pipeline (TWAS-Z + Mendelian Randomization)

Per-region pipeline that pairs QTL-derived weight vectors (TwasWeights and/or a QTL QtlFineMappingResult) with one or more GWAS studies (GwasSumStats) to produce per-tuple TWAS Z-scores and, when fine-mapping is supplied, Wald-ratio Mendelian Randomization estimates over the QTL credible sets.

Input combinations:

twasWeights alone (no fineMappingResult): TWAS Z only, no MR.
fineMappingResult alone (no twasWeights): TWAS Z derived from the SuSiE-style coefficients in each row's top loci (topLoci); plus MR.
both: TWAS Z computed from twasWeights; MR computed from fineMappingResult.

Usage

causalInferencePipeline(
  gwasSumStats,
  twasWeights = NULL,
  fineMappingResult = NULL,
  weightSelectionParam = WeightSelectionParam(),
  mrParam = MrParam(),
  combineMethods = NULL,
  alleleFlip = TRUE
)

Arguments

gwasSumStats: A GwasSumStats object. Must be QC'd (length(qcInfo(x)) > 0L).
twasWeights: Optional TwasWeights carrying per-(study, context, trait, method) weights. When supplied, drives the TWAS-Z computation.
fineMappingResult: Optional QtlFineMappingResult. When supplied, drives the MR computation and (when twasWeights = NULL) the TWAS-Z weights via the SuSiE-style coefficients on each entry's topLoci.
weightSelectionParam: How one TWAS method's weights are picked per (study, context, trait, gwasStudy) tuple, built with WeightSelectionParam: cutoff (minimum cross-validated r-squared, default 0), pvalCutoff (the CV-p-value gate, default Inf), and metric / pvalMetric naming which cvResult metrics those two read. Selection runs when cutoff > 0 or pvalCutoff is finite; it keeps the highest-metric method among those clearing both gates that also produced a finite TWAS Z. Otherwise every method is carried through. Ports the legacy twas_pipeline pick_best_model / update_twas_method.
mrParam: How Mendelian randomization is run, built with MrParam: method ("ivwPerVariant", the default, or "csAware"), the PIP threshold that method reads (pipCutoff or cpipCutoff), and pvalCutoff, the TWAS-p-value gate below which MR runs at all (default 1, no gate; otherwise the MR columns are NA).
combineMethods: Optional character vector forwarded to combinePValues for cross-method combination per (qtlStudy, context, trait, gwasStudy) group. NULL (default) skips combination.
alleleFlip: Logical, default TRUE. When TRUE, match QTL variants to the GWAS by (chrom, pos) with ref/alt swaps recognized and the exposure effect / weight sign-flipped accordingly; when FALSE, match on exact alleles only, so a ref/alt swap is treated as a distinct variant.

Value

A GRanges as described above.

LD reference compatibility check

If a QTL input (TwasWeights or QtlFineMappingResult) carries a non-NULL ldReference, it must be the same LD reference as the ldReference on gwasSumStats. 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: QC-ing the two sides separately trims each reference to its own surviving variants, and a partial overlap only warns. A mismatch of any of those, or no shared variant at all, is a hard error. A QTL input with ldReference = NULL (the fit was learned from individual-level data) skips the validation for that input.

Output shape

A long-format GRanges with one row per (qtlStudy, context, trait, method, gwasStudy) tuple, positioned at the variant span of the QTL weight set, with mcols:

qtlStudy, context, trait, method, gwasStudy: Identity columns.
twasZ, twasPval: Per-tuple TWAS Z and p-value.
waldRatio, waldRatioSe, mrPval: Per-tuple IVW-aggregated Wald-ratio MR estimate, standard error, and p-value. Present only when MR was computed; NA otherwise.
nIV: Number of instrumental variables used in the MR aggregation.

Examples

data(qtlDatasetExample)
data(gwasSumStatsS4Example)
tw <- twasWeightsPipeline(qtlDatasetExample,
  methods = "lasso", cisWindow = 1e6)
#>   Fitting lasso ...
#>   Fitting lasso done in 0.3s
causalInferencePipeline(
  gwasSumStats = gwasSumStatsS4Example, twasWeights = tw)
#> GRanges object with 1 range and 14 metadata columns:
#>       seqnames            ranges strand |    qtlStudy     context        trait
#>          <Rle>         <IRanges>  <Rle> | <character> <character>  <character>
#>   [1]    chr22 14850625-15826991      * |      study1       brain ENSG_example
#>            method   gwasStudy     twasZ  twasPval waldRatio waldRatioSe
#>       <character> <character> <numeric> <numeric> <numeric>   <numeric>
#>   [1]       lasso      trait1  -0.18457  0.853566        NA          NA
#>          mrPval       nIV         Q        I2       nCs
#>       <numeric> <integer> <numeric> <numeric> <integer>
#>   [1]        NA      <NA>        NA        NA      <NA>
#>   -------
#>   seqinfo: 1 sequence from an unspecified genome; no seqlengths