Retain fSuSiE curves and bands while storing each joint fit once

.github/environment/pixi.toml

@@ -9,7 +9,7 @@ platforms = [
 
 [tasks]
 devtools_document = "cd $GITHUB_WORKSPACE; R -e 'devtools::document()'"
-devtools_test = "cd $GITHUB_WORKSPACE; R -e 'devtools::test()'"
+devtools_test = "cd $GITHUB_WORKSPACE; R -e 'devtools::test(stop_on_failure = TRUE)'"
 codecov = "cd $GITHUB_WORKSPACE; R -e 'covr::codecov(quiet = FALSE)'"
 rcmdcheck = "cd $GITHUB_WORKSPACE; R -e 'rcmdcheck::rcmdcheck()'"
 bioccheck_git_clone = "cd $GITHUB_WORKSPACE; R -e 'BiocCheck::BiocCheckGitClone()'"

R/AllGenerics.R

@@ -2196,8 +2196,8 @@ setGeneric("getCredibleSetSummary", function(x, ...) {
 #' the fitted effect curve and its credible band over the functional grid, one
 #' row per (credible set, grid position). Wraps the upstream fSuSiE band
 #' computation (fsusieR's internal \code{update_cal_credible_band.susiF} +
-#' \code{get_fitted_effect}). Requires an UNtrimmed fit (the wavelet slots are
-#' dropped by trimming); degrades to zero rows for non-fSuSiE or trimmed fits.
+#' \code{get_fitted_effect}). Uses saved bands in trimmed fits, or computes
+#' them from an untrimmed fit. Returns zero rows when functional fields are absent.
 #'
 #' @param x A \code{FineMappingRow} or \code{FineMappingResultBase}.
 #' @param ... Ignored.
@@ -2218,8 +2218,8 @@ setGeneric("fsusieCredibleBand", function(x, ...) {
 #' zero (the genomic footprint of each fSuSiE effect), from the upstream
 #' \code{fsusieR::affected_reg}. Returns a \code{GRanges} with the credible-set
 #' label, its purity, and the effect \code{direction}
-#' (\code{"pos"}/\code{"neg"}; upstream discards it). Requires an UNtrimmed fit;
-#' empty otherwise.
+#' (\code{"pos"}/\code{"neg"}; upstream discards it). Uses retained functional
+#' summaries or an untrimmed fit; empty when those fields are absent.
 #'
 #' @param x A \code{FineMappingRow} or \code{FineMappingResultBase}.
 #' @param ... Ignored.

R/fineMappingRow.R

@@ -179,10 +179,10 @@ fineMappingRow <- function(variantIds, susieFit, topLoci, cvResult = NULL) {
 .isFsusieFit <- function(fit) {
     !is.null(fit) &&
         is.list(fit) &&
-        !is.null(fit$fitted_wc2) &&
         !is.null(fit$fitted_func) &&
         !is.null(fit$outing_grid) &&
-        !is.null(fit$alpha)
+        (!is.null(fit$cred_band) ||
+             (!is.null(fit$fitted_wc2) && !is.null(fit$alpha)))
 }
 
 # Populate `obj$cred_band` via fsusieR's wavethresh/GenW band computation. That
@@ -192,6 +192,9 @@ fineMappingRow <- function(variantIds, susieFit, topLoci, cvResult = NULL) {
 # an upstream fsusieR change surfaces as a clear error, not a silent NULL.
 # @noRd
 .fsusiePopulateCredibleBand <- function(fit) {
+    # Native post_processing="none" can leave zero-filled placeholder bands.
+    # Only reuse bands after trimming has removed the reconstruction fields.
+    if (is.null(fit$fitted_wc2) && !is.null(fit$cred_band)) return(fit)
     fn <- tryCatch(
         get("update_cal_credible_band.susiF", envir = asNamespace("fsusieR")),
         error = function(e) NULL

R/fineMappingWrappers.R

@@ -1857,11 +1857,15 @@ trimFinemappingFit <- function(fit, effectIdx, method, csTables) {
     if (method == "mvsusie") {
         trimmed <- .trimAddMvsusie(trimmed, fit, effectIdx)
     }
-    # fSuSiE: keep the precomputed variants x features TWAS weight matrix
-    # (fsusieWeights output, attached as $coef before trimming) so downstream
-    # TWAS can read it without the dropped wavelet slots.
-    if (method == "fsusie" && !is.null(fit$coef)) {
-        trimmed$coef <- fit$coef
+    # fSuSiE: retain functional summaries and precomputed weights. The base
+    # trimmed fit already retains alpha, lbf_variable and sets for coloc.
+    if (method == "fsusie") {
+        if (.isFsusieFit(fit)) {
+            fit <- .fsusiePopulateCredibleBand(fit)
+        }
+        keep <- intersect(c("coef", "fitted_func", "cred_band", "outing_grid",
+                            "cs", "csd_X", "lfsr_func"), names(fit))
+        trimmed[keep] <- fit[keep]
     }
     class(trimmed) <- unique(c(method, "susie"))
     trimmed

R/jointEngine.R

@@ -257,7 +257,7 @@ setMethod(
 )
 
 # fsusie joint fit (functional SuSiE over the trait domain; individual-level,
-# cross-trait). One per-condition entry per trait, with an optional CV slice.
+# cross-trait). Keep the first post-processed entry for this joint fit.
 # @noRd
 .jointFitFsusie <- function(group, Xc, Yc, nCond, cfg, args) {
     if (length(.jgTraitPos(group)) != nCond) {
@@ -275,22 +275,24 @@ setMethod(
     )
     fit <- exec(fitFsusie, !!!fitArgs)
     # Collapse the functional fit to a variants x features weight matrix now
-    # (trimming later drops fitted_wc/csd_X); store on $coef so a trimmed fit
+    # (trimming later drops fitted_wc); store on $coef so a trimmed fit
     # can still yield TWAS weights.
     fit$coef <- tryCatch(
         fsusieWeights(fsusieFit = fit, variantIds = colnames(Xc)),
         error = function(e) NULL
     )
     fit <- .setFinemappingFitClass(fit, "fsusie")
     cvM <- .jointFsusieCv(Xc, Yc, group, cfg, args, verbose)
-    map(
-        seq_len(nCond),
-        .jointFsusieEntry,
+    entry <- .jointFsusieEntry(
+        r = 1L,
         fit = fit,
         cvM = cvM,
         Xc = Xc,
         cfg = cfg
     )
+    entry@susieFit$trait_names <- as.character(.jgConditions(group)$trait)
+    entry@susieFit$trait_positions <- .jgTraitPos(group)
+    list(entry)
 }
 
 # Per-fold fsusie CV slice, or NULL when CV is disabled.

tests/testthat/test_fineMappingRow.R

@@ -967,7 +967,7 @@ test_that("getCredibleSetSummary aggregates across a collection with entry ident
 # test_fsusieAccessors.R when R/fsusieAccessors.R was folded in)
 # ===========================================================================
 
-.fsa_makeFit <- function() {
+.fsa_makeFit <- function(post_processing = "none") {
     set.seed(1)
     n <- 150L
     p <- 24L
@@ -990,7 +990,7 @@ test_that("getCredibleSetSummary aggregates across a collection with entry ident
         Y = Y,
         pos = seq_len(J),
         L = 5,
-        post_processing = "none",
+        post_processing = post_processing,
         verbose = FALSE
     ))
 }
@@ -1033,6 +1033,36 @@ test_that("fsusieCredibleBand returns a long effect + band table (lower <= effec
     expect_true(all(grepl("^fsusie_", cb$cs)))
 })
 
+test_that("trimmed fSuSiE survives serialization with curves, bands and coloc inputs", {
+    skip_if_not_installed("fsusieR")
+    skip_if_not_installed("wavethresh")
+    for (mode in c("none", "TI")) {
+        fit <- .fsa_makeFit(mode)
+        expectedFit <- pecotmr:::.fsusiePopulateCredibleBand(fit)
+        tables <- list(list(sets = list(cs = list())))
+        trimmed <- pecotmr:::trimFinemappingFit(
+            fit, seq_along(fit$alpha), "fsusie", tables
+        )
+        entry <- unserialize(serialize(.fsa_entry(trimmed), NULL))
+        saved <- entry@susieFit
+        expect_null(saved$fitted_wc)
+        expect_null(saved$fitted_wc2)
+        expect_identical(saved$fitted_func, fit$fitted_func)
+        expect_identical(saved$cred_band, expectedFit$cred_band)
+        expect_equal(.fmrRowFsusieCredibleBand(entry), .fmrRowFsusieCredibleBand(.fsa_entry(fit)))
+        expect_equal(.fmrRowFsusieAffectedRegions(entry), .fmrRowFsusieAffectedRegions(.fsa_entry(fit)))
+        lbf <- pecotmr:::.colocExtractLbfFromEntry(
+            entry, FALSE, NULL, 0.5, 1e-9
+        )
+        expected <- pecotmr:::.asLbfMatrix(fit)
+        colnames(expected) <- getVariantIds(entry)
+        expect_equal(lbf$lbf, expected)
+        colocResult <- coloc::coloc.bf_bf(lbf$lbf, lbf$lbf)
+        expect_true("PP.H4.abf" %in% names(colocResult$summary))
+    }
+})
+
+
 test_that("fsusieAffectedRegions returns a GRanges with cs / purity / direction", {
     skip_if_not_installed("fsusieR")
     skip_if_not_installed("wavethresh")

tests/testthat/test_jointEngine.R

@@ -845,7 +845,7 @@ test_that(".runJointCell: cross-trait twas -> per-trait weight vectors", {
     expect_false(is.matrix(getWeights(pecotmr:::.collectionEntry(res, 1L))))
 })
 
-test_that("fitJointGroup(Individual, Fm): fsusie returns one entry per trait", {
+test_that("fitJointGroup(Individual, Fm): fsusie keeps one entry for all traits", {
     set.seed(12)
     n <- 10L
     X <- matrix(
@@ -888,9 +888,12 @@ test_that("fitJointGroup(Individual, Fm): fsusie returns one entry per trait", {
 
     entries <- pecotmr:::fitJointGroup(grp, pipe, "fsusie", list())
     expect_type(entries, "list")
-    expect_length(entries, 2L) # one per trait
+    expect_length(entries, 1L)
     expect_s4_class(entries[[1L]], "FineMappingRow")
+    expect_identical(captured$Y, Y)
     expect_equal(captured$pos, c(100, 200)) # functional domain threaded
+    expect_identical(entries[[1L]]@susieFit$trait_names, c("G1", "G2"))
+    expect_equal(entries[[1L]]@susieFit$trait_positions, c(100, 200))
 })
 
 test_that("fitJointGroup(Individual, Fm): fsusie without pos errors; unknown token errors", {
@@ -1787,7 +1790,7 @@ test_that("fitJointGroup(Individual, Fm): fsusie honest per-fold CV is attached"
     )
     entries <- pecotmr:::fitJointGroup(grp, pipe, "fsusie", list())
     expect_true(cvCalled) # CV path exercised
-    expect_length(entries, 2L)
+    expect_length(entries, 1L)
 })
 
 test_that("fitJointGroup(Individual, Fm): SER pre-screen skips when < 2 survivors", {