diff --git a/DESCRIPTION b/DESCRIPTION index 86506d96..d5da4071 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -35,7 +35,7 @@ Imports: readr, rlang, stringr, - susieR (>= 0.16.2), + susieR (>= 0.16.4), tibble, tictoc, tidyr, diff --git a/R/fineMappingPipeline.R b/R/fineMappingPipeline.R index 50592e0b..9d37cfba 100644 --- a/R/fineMappingPipeline.R +++ b/R/fineMappingPipeline.R @@ -281,14 +281,14 @@ #' \code{rMismatch != "none"}) and \code{rFinite} is \code{NULL}, in which #' case it defaults to the LD-panel sample size \code{getNSamples(ldSketch)}. #' @param rMismatch \code{QtlSumStats} / \code{GwasSumStats} only. LD-mismatch -#' mode forwarded to \code{susieR::susie_rss()} as \code{R_mismatch} (e.g. -#' \code{"eb"} for empirical Bayes). Default \code{"none"} (susieR's default). -#' @param rMismatchMethod \code{QtlSumStats} / \code{GwasSumStats} only. Optional -#' \code{R_mismatch_method} forwarded to \code{susieR::susie_rss()} when -#' non-\code{NULL}. -#' @param checkPrior \code{QtlSumStats} / \code{GwasSumStats} only. Optional -#' \code{check_prior} forwarded to \code{susieR::susie_rss()} when -#' non-\code{NULL}. +#' correction mode forwarded to \code{susieR::susie_rss()} as +#' \code{R_mismatch}: \code{"none"} (default), \code{"eb"} (empirical Bayes), +#' or \code{"eb_mix"} (residual-mixture EB). +#' @param rssControl \code{QtlSumStats} / \code{GwasSumStats} only. Optional +#' named list of \code{susieR::susie_rss_control()} settings (e.g. +#' \code{check_prior}, \code{mismatch_estimator}), forwarded as +#' \code{susie_rss()}'s \code{control} argument. Default \code{NULL} leaves +#' the \code{susie_rss_control()} defaults in place. #' @param keepFullFit \code{QtlSumStats} / \code{GwasSumStats} only. Controls #' retention of the pre-fallback multi-effect SuSiE-RSS fit when #' \code{serFallback=TRUE}: \code{"fallback"} (default) keeps it only for @@ -1751,8 +1751,7 @@ setMethod("fineMappingPipeline", "QtlSumStats", serFallback = FALSE, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, - checkPrior = NULL, + rssControl = NULL, keepFullFit = "fallback", ...) { .fmAssertQcd(data) @@ -1872,8 +1871,8 @@ setMethod("fineMappingPipeline", "QtlSumStats", af = afByVar, fullFit = fullFit, fullFitAlphaOnly = fullFitAlphaOnly, includeAllCs = includeAllCs, serFallback = serFallback, rFinite = rFiniteResolved, - rMismatch = rMismatch, rMismatchMethod = rMismatchMethod, - checkPrior = checkPrior, keepFullFit = keepFullFit) + rMismatch = rMismatch, rssControl = rssControl, + keepFullFit = keepFullFit) # The method column carries the bare token, independent of the # postprocess class. for (tk in names(ents)) pushRow(st, ctx, tr, tk, ents[[tk]]) @@ -1953,8 +1952,7 @@ setMethod("fineMappingPipeline", "GwasSumStats", serFallback = FALSE, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, - checkPrior = NULL, + rssControl = NULL, keepFullFit = "fallback", ...) { .fmAssertQcd(data) @@ -2038,8 +2036,8 @@ setMethod("fineMappingPipeline", "GwasSumStats", af = afByVar, fullFit = fullFit, fullFitAlphaOnly = fullFitAlphaOnly, includeAllCs = includeAllCs, serFallback = serFallback, rFinite = rFiniteResolved, - rMismatch = rMismatch, rMismatchMethod = rMismatchMethod, - checkPrior = checkPrior, keepFullFit = keepFullFit) + rMismatch = rMismatch, rssControl = rssControl, + keepFullFit = keepFullFit) for (tk in names(ents)) pushRow(st, tk, region_id, ents[[tk]]) } diff --git a/R/fineMappingWrappers.R b/R/fineMappingWrappers.R index c66772c9..c663fe79 100644 --- a/R/fineMappingWrappers.R +++ b/R/fineMappingWrappers.R @@ -1992,7 +1992,7 @@ mergeSusieCs <- function(fineMappingResult, coverage = 0.95) { .fmFitSusieRss <- function(z, R, n, token, chainFromInf = NULL, coverage = 0.95, userArgs = NULL, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, checkPrior = NULL) { + rssControl = NULL) { info <- .fineMappingMethodCapabilities[[token]] if (is.null(info) || identical(info$unmappableEffects, NA_character_)) { stop(".fmFitSusieRss: token '", token, "' is not a SuSiE-family method.") @@ -2006,8 +2006,16 @@ mergeSusieCs <- function(fineMappingResult, coverage = 0.95) { # methodArgs (folded in after) still override them. baseArgs$R_finite <- rFinite baseArgs$R_mismatch <- rMismatch - if (!is.null(rMismatchMethod)) baseArgs$R_mismatch_method <- rMismatchMethod - if (!is.null(checkPrior)) baseArgs$check_prior <- checkPrior + # Optional susie_rss_control() settings (e.g. check_prior, mismatch_estimator), + # supplied as a named list and forwarded as susie_rss()'s `control` argument. + if (!is.null(rssControl)) { + if (!is.list(rssControl) || is.null(names(rssControl)) || + any(!nzchar(names(rssControl)))) { + stop(".fmFitSusieRss: `rssControl` must be a named list of ", + "susieR::susie_rss_control() settings.") + } + baseArgs$control <- do.call(susieR::susie_rss_control, rssControl) + } if (!is.null(chainFromInf) && token != "susieInf") { # SuSiE-RSS(-ash) initialised from a SuSiE-inf fit; userArgs folded into the # arg prep so L_greedy is clamped rather than passed through raw. @@ -2115,8 +2123,8 @@ mergeSusieCs <- function(fineMappingResult, coverage = 0.95) { fullFit = FALSE, fullFitAlphaOnly = TRUE, includeAllCs = FALSE, serFallback = FALSE, rFinite = NULL, - rMismatch = "none", rMismatchMethod = NULL, - checkPrior = NULL, keepFullFit = "fallback") { + rMismatch = "none", rssControl = NULL, + keepFullFit = "fallback") { chainLocal <- .fmResolveSusieChain(toRun, addSusieInf) infFit <- NULL if (chainLocal$runInf) { @@ -2125,8 +2133,7 @@ mergeSusieCs <- function(fineMappingResult, coverage = 0.95) { infFit <- .fmFitSusieRss(z, R, n, "susieInf", coverage = coverage, userArgs = methodArgs[["susieInf"]], rFinite = rFinite, rMismatch = rMismatch, - rMismatchMethod = rMismatchMethod, - checkPrior = checkPrior) + rssControl = rssControl) } out <- list() for (tk in toRun) { @@ -2156,8 +2163,7 @@ mergeSusieCs <- function(fineMappingResult, coverage = 0.95) { fit <- .fmFitSusieRss(z, R, n, tk, chainFromInf = chainFrom, coverage = coverage, userArgs = methodArgs[[tk]], rFinite = rFinite, rMismatch = rMismatch, - rMismatchMethod = rMismatchMethod, - checkPrior = checkPrior) + rssControl = rssControl) rfd <- fit$R_finite_diagnostics flag <- if (!is.null(rfd) && !is.null(rfd$R_reliability_flag)) isTRUE(rfd$R_reliability_flag) else NA diff --git a/man/fineMappingPipeline.Rd b/man/fineMappingPipeline.Rd index 4508561a..362c0d0b 100644 --- a/man/fineMappingPipeline.Rd +++ b/man/fineMappingPipeline.Rd @@ -117,8 +117,7 @@ fineMappingPipeline(data, ...) serFallback = FALSE, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, - checkPrior = NULL, + rssControl = NULL, keepFullFit = "fallback", ... ) @@ -143,8 +142,7 @@ fineMappingPipeline(data, ...) serFallback = FALSE, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, - checkPrior = NULL, + rssControl = NULL, keepFullFit = "fallback", ... ) @@ -355,16 +353,15 @@ except when a finite/EB mode is active (\code{serFallback=TRUE} or case it defaults to the LD-panel sample size \code{getNSamples(ldSketch)}.} \item{rMismatch}{\code{QtlSumStats} / \code{GwasSumStats} only. LD-mismatch -mode forwarded to \code{susieR::susie_rss()} as \code{R_mismatch} (e.g. -\code{"eb"} for empirical Bayes). Default \code{"none"} (susieR's default).} - -\item{rMismatchMethod}{\code{QtlSumStats} / \code{GwasSumStats} only. Optional -\code{R_mismatch_method} forwarded to \code{susieR::susie_rss()} when -non-\code{NULL}.} - -\item{checkPrior}{\code{QtlSumStats} / \code{GwasSumStats} only. Optional -\code{check_prior} forwarded to \code{susieR::susie_rss()} when -non-\code{NULL}.} +correction mode forwarded to \code{susieR::susie_rss()} as +\code{R_mismatch}: \code{"none"} (default), \code{"eb"} (empirical Bayes), +or \code{"eb_mix"} (residual-mixture EB).} + +\item{rssControl}{\code{QtlSumStats} / \code{GwasSumStats} only. Optional +named list of \code{susieR::susie_rss_control()} settings (e.g. +\code{check_prior}, \code{mismatch_estimator}), forwarded as +\code{susie_rss()}'s \code{control} argument. Default \code{NULL} leaves +the \code{susie_rss_control()} defaults in place.} \item{keepFullFit}{\code{QtlSumStats} / \code{GwasSumStats} only. Controls retention of the pre-fallback multi-effect SuSiE-RSS fit when diff --git a/tests/testthat/test_fineMappingPipeline.R b/tests/testthat/test_fineMappingPipeline.R index ea011b3c..857cef68 100644 --- a/tests/testthat/test_fineMappingPipeline.R +++ b/tests/testthat/test_fineMappingPipeline.R @@ -131,7 +131,7 @@ context("fineMappingPipeline") .fmp_mockFitRss <- function() { function(z, R, n, token, chainFromInf = NULL, coverage = 0.95, userArgs = NULL, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, checkPrior = NULL) { + rssControl = NULL) { list(token = token, n_variants = length(z)) } } @@ -143,12 +143,11 @@ context("fineMappingPipeline") .fmp_mockFitRssDiag <- function(flag = FALSE, capture = NULL) { function(z, R, n, token, chainFromInf = NULL, coverage = 0.95, userArgs = NULL, rFinite = NULL, rMismatch = "none", - rMismatchMethod = NULL, checkPrior = NULL) { + rssControl = NULL) { if (!is.null(capture)) { capture$rFinite <- rFinite capture$rMismatch <- rMismatch - capture$rMismatchMethod <- rMismatchMethod - capture$checkPrior <- checkPrior + capture$rssControl <- rssControl } fit <- list(token = token, n_variants = length(z), multiTag = "multi") fit$R_finite_diagnostics <- list( @@ -1602,6 +1601,42 @@ test_that("fineMappingPipeline(GwasSumStats): rFinite/rMismatch forwarded; rFini fineMappingPipeline(gss, methods = "susie", addSusieInf = FALSE, serFallback = TRUE, rFinite = 12345, rMismatch = "eb")) expect_equal(cap2$rFinite, 12345) + + # eb_mix (susieR >= 0.16.6) is passed through to R_mismatch unrestricted. + cap3 <- new.env(parent = emptyenv()) + local_mocked_bindings( + extractBlockGenotypes = .fmp_mockExtractor(), + .fmFitSusieRss = .fmp_mockFitRssDiag(flag = FALSE, capture = cap3), + .fmPostprocessOne = .fmp_mockPostprocess(), + .package = "pecotmr") + suppressMessages( + fineMappingPipeline(gss, methods = "susie", addSusieInf = FALSE, + serFallback = TRUE, rMismatch = "eb_mix")) + expect_equal(cap3$rMismatch, "eb_mix") + + # rssControl (susie_rss_control() settings) forwarded through to the fitter. + cap4 <- new.env(parent = emptyenv()) + local_mocked_bindings( + extractBlockGenotypes = .fmp_mockExtractor(), + .fmFitSusieRss = .fmp_mockFitRssDiag(flag = FALSE, capture = cap4), + .fmPostprocessOne = .fmp_mockPostprocess(), + .package = "pecotmr") + suppressMessages( + fineMappingPipeline(gss, methods = "susie", addSusieInf = FALSE, + rssControl = list(check_prior = TRUE, + mismatch_estimator = "map"))) + expect_equal(cap4$rssControl, list(check_prior = TRUE, + mismatch_estimator = "map")) + # default: no rssControl -> fitter sees NULL + cap5 <- new.env(parent = emptyenv()) + local_mocked_bindings( + extractBlockGenotypes = .fmp_mockExtractor(), + .fmFitSusieRss = .fmp_mockFitRssDiag(flag = FALSE, capture = cap5), + .fmPostprocessOne = .fmp_mockPostprocess(), + .package = "pecotmr") + suppressMessages( + fineMappingPipeline(gss, methods = "susie", addSusieInf = FALSE)) + expect_null(cap5$rssControl) }) test_that("fineMappingPipeline(GwasSumStats): keepFullFit='all' retains fit on non-fallback region", { diff --git a/tests/testthat/test_fineMappingWrappers.R b/tests/testthat/test_fineMappingWrappers.R index 6431fec6..278e1200 100644 --- a/tests/testthat/test_fineMappingWrappers.R +++ b/tests/testthat/test_fineMappingWrappers.R @@ -1768,6 +1768,36 @@ test_that("fitMvsusie forwards arguments to mvsusieR::mvsusie", { expect_equal(r$coverage, 0.9) }) +test_that(".fmFitSusieRss rejects a non-named-list rssControl", { + z <- c(0.5, -1.2, 2.0); R <- diag(3) + expect_error( + pecotmr:::.fmFitSusieRss(z, R, n = 1000, token = "susie", rssControl = "nope"), + "named list", fixed = TRUE) + expect_error( + pecotmr:::.fmFitSusieRss(z, R, n = 1000, token = "susie", rssControl = list(1, 2)), + "named list", fixed = TRUE) + expect_error( # partially named + pecotmr:::.fmFitSusieRss(z, R, n = 1000, token = "susie", + rssControl = list(check_prior = TRUE, 2)), + "named list", fixed = TRUE) +}) + +test_that(".fmFitSusieRss forwards rssControl to susie_rss as control", { + cap <- new.env(parent = emptyenv()) + local_mocked_bindings( + susie_rss_control = function(...) list(.tag = "ctrl", ...), + susie_rss = function(...) { args <- list(...); cap$control <- args$control + list(sets = list(), pip = numeric()) }, + .package = "susieR") + z <- c(0.5, -1.2, 2.0); R <- diag(3) + pecotmr:::.fmFitSusieRss(z, R, n = 1000, token = "susie", + rssControl = list(check_prior = TRUE, + mismatch_estimator = "map")) + expect_equal(cap$control$.tag, "ctrl") + expect_true(cap$control$check_prior) + expect_equal(cap$control$mismatch_estimator, "map") +}) + test_that("fitMvsusieRss forwards arguments to mvsusieR::mvsusie_rss", { skip_if_not_installed("mvsusieR") local_mocked_bindings(