Hierarchical Multiple-Testing Correction for cis-QTL Association

Apply hierarchical (local per-gene + global across-gene) multiple-testing correction to a per-gene QtlSumStats of cis-QTL association statistics (one row per gene/trait), returning the same object enriched with the corrected-statistic columns. See qtlAssociationPostprocess.

Usage

qtlAssociationPostprocess(x, ...)

# S4 method for class 'QtlSumStats'
qtlAssociationPostprocess(
  x,
  fdrThreshold = 0.05,
  mafCutoff = 0,
  cisWindow = 0,
  methods = c("permutation", "bonferroni"),
  pvalueCol = "P",
  afCol = "af"
)

Arguments

x

A QtlSumStats.

...

Correction arguments.

fdrThreshold

Event- and variant-level FDR threshold (default 0.05).

mafCutoff

Minor-allele-frequency cutoff for the FILTERED Bonferroni flavour (fold of the per-variant af mcol). 0 disables it.

cisWindow

cis-window (bp) for the FILTERED Bonferroni flavour, applied to the per-variant tss_distance/tes_distance mcols. 0 disables it. When either mafCutoff or cisWindow is > 0 the *_filtered columns are produced (requires the n_variants_filtered row column).

methods

Correction families to compute: any of "permutation" (needs p_beta/beta_shape1/beta_shape2) and "bonferroni" (needs n_variants). The q-value SNP method adds no stored column – it is a significance query (see getSignificantQtls).

pvalueCol, afCol

Entry mcol names for the per-variant p-value / allele frequency (defaults "P" / "af").

Value

The input QtlSumStats with added correction columns.

Examples

pBeta <- c(10^(-c(8, 7, 6, 5, 4, 3)), stats::ppoints(54))
G <- length(pBeta)
entries <- lapply(seq_len(G), function(i) {
  gr <- GenomicRanges::GRanges("chr1",
    IRanges::IRanges(seq(1000L, by = 50L, length.out = 4L), width = 1L))
  S4Vectors::mcols(gr) <- S4Vectors::DataFrame(
    SNP = paste0("g", i, "_v", 1:4), A1 = "A", A2 = "G",
    P = c(pBeta[i] / 5, 0.2, 0.5, 0.8),
    af = c(0.3, 0.2, 0.005, 0.4),
    tss_distance = c(0L, 500L, 900000L, 2000000L),
    tes_distance = c(0L, 500L, 900000L, 2000000L))
  gr
})
qss <- QtlSumStats(study = rep("s", G), context = rep("brain", G),
  trait = paste0("g", seq_len(G)), entry = entries, genome = "hg19",
  n_variants = rep(50L, G), n_variants_filtered = rep(30L, G),
  p_beta = pBeta, beta_shape1 = rep(1, G), beta_shape2 = rep(200, G))
qtlAssociationPostprocess(qss)
#> QtlSumStats: 60 entries, genome build hg19
#>   1 studies, 1 contexts, 60 traits
#>   LD sketch: none (LD-free)