Unify compact top-loci export and SuSiE-Inf as the option for univaraite pipeline

NAMESPACE

@@ -25,6 +25,7 @@ export(bayes_c_weights)
 export(bayes_l_weights)
 export(bayes_n_weights)
 export(bayes_r_weights)
+export(build_top_loci)
 export(check_ld)
 export(clean_context_names)
 export(coloc_post_processor)

R/multivariate_pipeline.R

@@ -71,6 +71,7 @@ multivariate_analysis_pipeline <- function(
     maf,
     X_variance = NULL,
     other_quantities = list(),
+    region = NULL,
     # filters
     imiss_cutoff = 1.0,
     maf_cutoff = 0.01,
@@ -304,7 +305,8 @@ multivariate_analysis_pipeline <- function(
     secondary_coverage = sec_coverage,
     signal_cutoff = signal_cutoff,
     min_abs_corr = min_abs_corr,
-    other_quantities = other_quantities
+    other_quantities = other_quantities,
+    region = region
   )
   res <- c(res, format_finemapping_output(mvsusie_post, primary_method = "mvsusie"))
   res$total_time_elapsed <- proc.time() - st

R/univariate_pipeline.R

@@ -1,8 +1,8 @@
 #' Univariate Analysis Pipeline
 #'
 #' This function performs univariate analysis for fine-mapping and Transcriptome-Wide Association Study (TWAS)
-#' with optional cross-validation. Fine-mapping fits SuSiE-inf first and then
-#' fits SuSiE initialized from the SuSiE-inf result.
+#' with optional cross-validation. By default, fine-mapping fits SuSiE-inf first
+#' and then fits SuSiE initialized from the SuSiE-inf result.
 #'
 #' @param X A matrix of genotype data where rows represent samples and columns represent genetic variants.
 #' @param Y A vector of phenotype measurements.
@@ -11,6 +11,7 @@
 #' @param maf A vector of minor allele frequencies for each variant in X.
 #' @param X_variance Optional variance of X. Default is NULL.
 #' @param other_quantities A list of other quantities to be carried into fine-mapping post-processing. Default is an empty list.
+#' @param region Optional \code{"chr:start-end"} string for the analysis region. Default is NULL.
 #' @param imiss_cutoff Individual missingness cutoff. Default is 1.0.
 #' @param maf_cutoff Minor allele frequency cutoff. Default is NULL.
 #' @param xvar_cutoff Variance cutoff for X. Default is 0.05.
@@ -25,6 +26,9 @@
 #' @param finemapping_extra_opts Additional options passed to \code{susieR::susie()}.
 #'   SuSiE-inf is always fitted with \code{refine = FALSE}; the ordinary SuSiE
 #'   fit keeps these options and is initialized with \code{model_init}.
+#' @param add_susie_inf Whether to fit SuSiE-inf before ordinary SuSiE. Default
+#'   is TRUE for existing pipeline compatibility. If FALSE, ordinary SuSiE is
+#'   fitted directly and SuSiE-inf fitted objects/results are not returned.
 #' @param twas_weights Whether to compute TWAS weights. Default is TRUE.
 #' @param sample_partition Optional data frame with Sample and Fold columns for cross-validation. Default is NULL.
 #' @param max_cv_variants The maximum number of variants to be included in cross-validation. Default is -1 (no limit).
@@ -44,6 +48,7 @@ univariate_analysis_pipeline <- function(
     Y_scalar = 1,
     X_variance = NULL,
     other_quantities = list(),
+    region = NULL,
     # filters
     imiss_cutoff = 1.0,
     maf_cutoff = NULL,
@@ -58,6 +63,7 @@ univariate_analysis_pipeline <- function(
     coverage = c(0.95, 0.7, 0.5),
     min_abs_corr = 0.8,
     finemapping_extra_opts = list(refine = TRUE),
+    add_susie_inf = TRUE,
     # TWAS weights and CV for TWAS weights
     twas_weights = TRUE,
     sample_partition = NULL,
@@ -75,6 +81,12 @@ univariate_analysis_pipeline <- function(
   if (!is.numeric(Y_scalar) || length(Y_scalar) != 1) stop("Y_scalar must be a numeric scalar")
   if (!is.numeric(L) || L <= 0) stop("L must be a positive integer")
   if (!is.null(L_greedy) && (!is.numeric(L_greedy) || L_greedy <= 0)) stop("L_greedy must be NULL or a positive integer")
+  if (!is.logical(add_susie_inf) || length(add_susie_inf) != 1 || is.na(add_susie_inf)) {
+    stop("add_susie_inf must be TRUE or FALSE")
+  }
+  if (!isTRUE(add_susie_inf) && isTRUE(twas_weights)) {
+    stop("add_susie_inf = FALSE is not compatible with twas_weights = TRUE")
+  }
 
   # Initial PIP check
   if (pip_cutoff_to_skip != 0) {
@@ -112,17 +124,28 @@ univariate_analysis_pipeline <- function(
   st <- proc.time()
   res <- list()
 
-  message("Fitting SuSiE-inf model on input data ...")
-  message("Fitting SuSiE model initialized by SuSiE-inf ...")
-  fitted_models <- fit_susie_inf_then_susie(
-    X,
-    Y,
-    args = modifyList(
-      finemapping_extra_opts,
-      list(L = L, L_greedy = L_greedy, coverage = coverage[1])
-    )
+  susie_args <- modifyList(
+    finemapping_extra_opts,
+    list(L = L, L_greedy = L_greedy, coverage = coverage[1])
   )
-  res$susie_inf_fitted <- fitted_models[["susie_inf"]]
+  if (isTRUE(add_susie_inf)) {
+    message("Fitting SuSiE-inf model on input data ...")
+    message("Fitting SuSiE model initialized by SuSiE-inf ...")
+    fitted_models <- fit_susie_inf_then_susie(
+      X,
+      Y,
+      args = susie_args
+    )
+    res$susie_inf_fitted <- fitted_models[["susie_inf"]]
+  } else {
+    message("Fitting SuSiE model on input data ...")
+    fitted_models <- list(
+      susie = .set_finemapping_fit_class(
+        do.call(susie, c(list(X = X, y = Y), susie_args)),
+        "susie"
+      )
+    )
+  }
   res$susie_fitted <- fitted_models[["susie"]]
 
   # Process SuSiE results
@@ -137,10 +160,13 @@ univariate_analysis_pipeline <- function(
     secondary_coverage = if (length(coverage) > 1) coverage[-1] else NULL,
     signal_cutoff = signal_cutoff,
     min_abs_corr = min_abs_corr,
-    other_quantities = other_quantities
+    other_quantities = other_quantities,
+    region = region
   )
   res <- c(res, format_finemapping_output(susie_post, primary_method = "susie"))
-  res$susie_inf_result_trimmed <- susie_post$finemapping_results$susie_inf$result_trimmed
+  if (!is.null(susie_post$finemapping_results$susie_inf)) {
+    res$susie_inf_result_trimmed <- susie_post$finemapping_results$susie_inf$result_trimmed
+  }
   res$total_time_elapsed <- proc.time() - st
 
   # TWAS weights and cross-validation