sLDSC Postprocessing Pipeline
Source:R/sldscPostprocessingPipeline.R
sldscPostprocessingPipeline.RdPostprocess polyfun's per-trait sLDSC outputs (already loaded
into an SldscData object) into a single results object with
per-trait tau*, EnrichStat with back-solved jackknife SE, and a
DerSimonian-Laird random-effects meta-analysis across traits. All file I/O
is done up front by the reader functions (readSldscAnnot,
readSldscFrq, readSldscTrait); this pipeline is
pure computation over the in-memory SldscData.
Usage
sldscPostprocessingPipeline(
sldscData,
mafCutoff = 0.05,
targetCategories = NULL,
targetLabels = NULL
)Arguments
- sldscData
An
SldscDataobject bundling the annotation table, the reference-panel allele frequencies, and the per-trait single/joint polyfun runs.- mafCutoff
Numeric MAF cutoff applied via the object's frq table. Default
0.05. Set to0to opt out (requires frq data when> 0).- targetCategories
Optional character vector of target annotation names to retain. Auto-detected from the joint run (or first single run) when
NULL.- targetLabels
Optional display names, same length / order as
targetCategories, applied to every output column / tau* block colname.
Value
A list with per_trait (per-trait standardised tables), meta
tables (tauStar, enrichment, enrichstat), and a
params record of the call options.
Examples
mkRun <- function(cats) {
n <- length(cats)
list(categories = cats, tau = setNames(rep(1e-7, n), cats),
tauSe = setNames(rep(3e-8, n), cats),
enrichment = setNames(rep(2, n), cats),
enrichmentSe = setNames(rep(0.4, n), cats),
enrichmentP = setNames(rep(0.01, n), cats),
propH2 = setNames(rep(0.2, n), cats),
propSnps = setNames(rep(0.1, n), cats), h2g = 0.3,
annotationJackknifeCoefs = matrix(1e-7 * seq(0.9, 1.1, length.out = 10),
10, n, dimnames = list(NULL, cats)),
nBlocks = 10L)
}
annot <- data.frame(CHR = c(1, 1, 1, 2, 2, 2), SNP = paste0("rs", 1:6),
annot_A = c(1, 0, 1, 0, 1, 0), annot_B = c(2.1, 1.8, 2.5, 1.9, 2.3, 2))
frq <- data.frame(CHR = c(1, 1, 1, 2, 2, 2), SNP = paste0("rs", 1:6),
MAF = rep(0.2, 6))
mkTrait <- function() {
list(single = list(mkRun(c("annot_A_0", "baselineLD_0")),
mkRun(c("annot_B_0", "baselineLD_0"))),
joint = mkRun(c("annot_A_0", "annot_B_0", "baselineLD_0")))
}
traits <- setNames(list(mkTrait(), mkTrait()), c("traitX", "traitY"))
sldscData <- SldscData(annotData = annot, frqData = frq, traits = traits)
sldscPostprocessingPipeline(sldscData = sldscData)
#> [sldsc] Computing M_ref...
#> [sldsc] M_ref = 6 (MAF cutoff 0.05)
#> [sldsc] Computing per-annotation sd...
#> [sldsc] sd computed for 2 annotation columns
#> [sldsc] Detecting binary vs continuous annotations...
#> [sldsc] Auto-detected 2 target categories
#> [sldsc] Detected 1 baseline annotations: baselineLD_0
#> [sldsc] Standardizing 2 traits...
#> [sldsc] Running random-effects meta across traits...
#> $per_trait
#> $per_trait$traitX
#> $per_trait$traitX$summary
#> # A tibble: 2 × 20
#> target isBinary tauSingle tauSeSingle tauStarSingle tauStarSeSingle
#> <chr> <lgl> <dbl> <dbl> <dbl> <dbl>
#> 1 annot_A_0 TRUE 0.0000001 0.00000003 0.00000115 0.000000221
#> 2 annot_B_0 FALSE 0.0000001 0.00000003 0.000000577 0.000000111
#> # ℹ 14 more variables: enrichmentSingle <dbl>, enrichmentSeSingle <dbl>,
#> # enrichmentPSingle <dbl>, enrichstatSingle <dbl>, enrichstatSeSingle <dbl>,
#> # tauJoint <dbl>, tauSeJoint <dbl>, tauStarJoint <dbl>, tauStarSeJoint <dbl>,
#> # enrichmentJoint <dbl>, enrichmentSeJoint <dbl>, enrichmentPJoint <dbl>,
#> # enrichstatJoint <dbl>, enrichstatSeJoint <dbl>
#>
#> $per_trait$traitX$tau_star_blocks_single
#> annot_A_0 annot_B_0
#> [1,] 1.039230e-06 5.196152e-07
#> [2,] 1.064890e-06 5.324452e-07
#> [3,] 1.090551e-06 5.452753e-07
#> [4,] 1.116211e-06 5.581053e-07
#> [5,] 1.141871e-06 5.709353e-07
#> [6,] 1.167531e-06 5.837653e-07
#> [7,] 1.193191e-06 5.965953e-07
#> [8,] 1.218851e-06 6.094253e-07
#> [9,] 1.244511e-06 6.222553e-07
#> [10,] 1.270171e-06 6.350853e-07
#>
#> $per_trait$traitX$tau_star_blocks_joint
#> annot_A_0 annot_B_0
#> [1,] 1.039230e-06 5.196152e-07
#> [2,] 1.064890e-06 5.324452e-07
#> [3,] 1.090551e-06 5.452753e-07
#> [4,] 1.116211e-06 5.581053e-07
#> [5,] 1.141871e-06 5.709353e-07
#> [6,] 1.167531e-06 5.837653e-07
#> [7,] 1.193191e-06 5.965953e-07
#> [8,] 1.218851e-06 6.094253e-07
#> [9,] 1.244511e-06 6.222553e-07
#> [10,] 1.270171e-06 6.350853e-07
#>
#> $per_trait$traitX$h2g
#> [1] 0.3
#>
#> $per_trait$traitX$nBlocks
#> [1] 10
#>
#>
#> $per_trait$traitY
#> $per_trait$traitY$summary
#> # A tibble: 2 × 20
#> target isBinary tauSingle tauSeSingle tauStarSingle tauStarSeSingle
#> <chr> <lgl> <dbl> <dbl> <dbl> <dbl>
#> 1 annot_A_0 TRUE 0.0000001 0.00000003 0.00000115 0.000000221
#> 2 annot_B_0 FALSE 0.0000001 0.00000003 0.000000577 0.000000111
#> # ℹ 14 more variables: enrichmentSingle <dbl>, enrichmentSeSingle <dbl>,
#> # enrichmentPSingle <dbl>, enrichstatSingle <dbl>, enrichstatSeSingle <dbl>,
#> # tauJoint <dbl>, tauSeJoint <dbl>, tauStarJoint <dbl>, tauStarSeJoint <dbl>,
#> # enrichmentJoint <dbl>, enrichmentSeJoint <dbl>, enrichmentPJoint <dbl>,
#> # enrichstatJoint <dbl>, enrichstatSeJoint <dbl>
#>
#> $per_trait$traitY$tau_star_blocks_single
#> annot_A_0 annot_B_0
#> [1,] 1.039230e-06 5.196152e-07
#> [2,] 1.064890e-06 5.324452e-07
#> [3,] 1.090551e-06 5.452753e-07
#> [4,] 1.116211e-06 5.581053e-07
#> [5,] 1.141871e-06 5.709353e-07
#> [6,] 1.167531e-06 5.837653e-07
#> [7,] 1.193191e-06 5.965953e-07
#> [8,] 1.218851e-06 6.094253e-07
#> [9,] 1.244511e-06 6.222553e-07
#> [10,] 1.270171e-06 6.350853e-07
#>
#> $per_trait$traitY$tau_star_blocks_joint
#> annot_A_0 annot_B_0
#> [1,] 1.039230e-06 5.196152e-07
#> [2,] 1.064890e-06 5.324452e-07
#> [3,] 1.090551e-06 5.452753e-07
#> [4,] 1.116211e-06 5.581053e-07
#> [5,] 1.141871e-06 5.709353e-07
#> [6,] 1.167531e-06 5.837653e-07
#> [7,] 1.193191e-06 5.965953e-07
#> [8,] 1.218851e-06 6.094253e-07
#> [9,] 1.244511e-06 6.222553e-07
#> [10,] 1.270171e-06 6.350853e-07
#>
#> $per_trait$traitY$h2g
#> [1] 0.3
#>
#> $per_trait$traitY$nBlocks
#> [1] 10
#>
#>
#>
#> $meta
#> $meta$tauStar
#> # A tibble: 2 × 9
#> target isBinary singleMean singleSe singleP nTraits jointMean jointSe
#> <chr> <lgl> <dbl> <dbl> <dbl> <int> <dbl> <dbl>
#> 1 annot_A_0 TRUE 0.00000115 0.000000156 1.52e-13 2 1.15e-6 1.56e-7
#> 2 annot_B_0 FALSE 0.000000577 0.0000000782 1.52e-13 2 5.77e-7 7.82e-8
#> # ℹ 1 more variable: jointP <dbl>
#>
#> $meta$enrichment
#> # A tibble: 2 × 6
#> target isBinary singleMean singleSe singleP nTraits
#> <chr> <lgl> <dbl> <dbl> <dbl> <int>
#> 1 annot_A_0 TRUE 2 0.283 0.000270 2
#> 2 annot_B_0 FALSE 2 0.283 0.000270 2
#>
#> $meta$enrichstat
#> # A tibble: 2 × 6
#> target isBinary singleMean singleSe singleP nTraits
#> <chr> <lgl> <dbl> <dbl> <dbl> <int>
#> 1 annot_A_0 TRUE 0.0556 0.0153 0.000270 2
#> 2 annot_B_0 FALSE 0.0556 0.0153 0.000270 2
#>
#>
#> $params
#> $params$maf_cutoff
#> [1] 0.05
#>
#> $params$M_ref
#> [1] 6
#>
#> $params$target_categories
#> [1] "annot_A_0" "annot_B_0"
#>
#> $params$n_baseline
#> [1] 1
#>
#> $params$baseline_categories
#> [1] "baselineLD_0"
#>
#> $params$trait_names
#> [1] "traitX" "traitY"
#>
#>