diff --git a/R/multivariate_pipeline.R b/R/multivariate_pipeline.R index 67f49715..079560d1 100644 --- a/R/multivariate_pipeline.R +++ b/R/multivariate_pipeline.R @@ -286,9 +286,9 @@ multivariate_analysis_pipeline <- function( message("Fitting mvSuSiE model on input data ...") res$mvsusie_fitted <- mvsusieR::mvsusie(X, Y, L = max_L, prior_variance = mvsusie_reweighted_mixture_prior$data_driven_prior_matrices, - residual_variance = resid_Y, precompute_covariances = FALSE, compute_objective = TRUE, - estimate_residual_variance = FALSE, estimate_prior_variance = TRUE, estimate_prior_method = "EM", - max_iter = mvsusie_max_iter, n_thread = 1, approximate = FALSE, verbosity = verbose, coverage = coverage[1] + residual_variance = resid_Y, estimate_residual_variance = FALSE, + max_iter = mvsusie_max_iter, + verbose = verbose, coverage = coverage[1] ) # Process mvSuSiE results diff --git a/R/regularized_regression.R b/R/regularized_regression.R index e63ff495..3e92a5ac 100644 --- a/R/regularized_regression.R +++ b/R/regularized_regression.R @@ -338,7 +338,9 @@ mrmash_weights <- function(mrmash_fit = NULL, X = NULL, Y = NULL, ...) { } #' @export -mvsusie_weights <- function(mvsusie_fit = NULL, X = NULL, Y = NULL, prior_variance = NULL, residual_variance = NULL, L = 30, ...) { +mvsusie_weights <- function(mvsusie_fit = NULL, X = NULL, Y = NULL, + prior_variance = NULL, residual_variance = NULL, + L = 30, verbose = FALSE, ...) { if (!requireNamespace("mvsusieR", quietly = TRUE)) { stop("Package 'mvsusieR' is required. Install with: devtools::install_github('stephenslab/mvsusieR')") } @@ -357,9 +359,9 @@ mvsusie_weights <- function(mvsusie_fit = NULL, X = NULL, Y = NULL, prior_varian mvsusie_fit <- mvsusieR::mvsusie( X = X, Y = Y, L = L, prior_variance = prior_variance, - residual_variance = residual_variance, precompute_covariances = FALSE, - compute_objective = TRUE, estimate_residual_variance = FALSE, estimate_prior_variance = TRUE, - estimate_prior_method = "EM", approximate = FALSE, ... + residual_variance = residual_variance, + estimate_residual_variance = FALSE, + verbose = verbose, ... ) } return(mvsusieR::coef.mvsusie(mvsusie_fit)[-1, ]) diff --git a/R/susie_wrapper.R b/R/susie_wrapper.R index b2ac41aa..a2144028 100644 --- a/R/susie_wrapper.R +++ b/R/susie_wrapper.R @@ -150,9 +150,8 @@ susie_wrapper <- function(X, y, init_L = 5, max_L = 30, l_step = 5, ...) { #' @param L Initial number of causal configurations. #' @param max_L Maximum number of causal configurations. #' @param l_step Step size for increasing L when the limit is reached. -#' @param stochastic_ld_sample Stochastic LD parameter passed to susie_rss. -#' NULL (default): no variance inflation. TRUE: infer sketch size from X -#' (requires X, not R). Integer: explicit sketch size (for R only). +#' @param R_finite Controls variance inflation to account for finite reference LD. +#' Passed to \code{susieR::susie_rss()}. #' @param ... Extra parameters passed to susie_rss (e.g., var_y, coverage). #' @return SuSiE RSS fit object after dynamic L adjustment #' @importFrom susieR susie_rss @@ -160,14 +159,14 @@ susie_wrapper <- function(X, y, init_L = 5, max_L = 30, l_step = 5, ...) { susie_rss_wrapper <- function(z, R = NULL, X = NULL, n = NULL, L = 10, max_L = 30, l_step = 5, coverage = 0.95, - stochastic_ld_sample = NULL, ...) { + R_finite = NULL, ...) { # Validate: exactly one of R or X if (is.null(R) && is.null(X)) stop("Either R or X must be provided.") if (!is.null(R) && !is.null(X)) stop("Only one of R or X should be provided, not both.") # Build argument list for susie_rss base_args <- list(z = z, n = n, L = L, coverage = coverage, - stochastic_ld_sample = stochastic_ld_sample, ...) + R_finite = R_finite, ...) if (!is.null(X)) base_args$X <- X else base_args$R <- R run_susie <- function(args) do.call(susie_rss, args) @@ -211,7 +210,8 @@ susie_rss_wrapper <- function(z, R = NULL, X = NULL, n = NULL, #' @param secondary_coverage Secondary coverage levels (default: c(0.7, 0.5)). #' @param signal_cutoff PIP cutoff for susie_post_processor (default: 0.1). #' @param min_abs_corr Minimum absolute correlation for CS purity (default: 0.8). -#' @param stochastic_ld_sample Passed to susie_rss. NULL, TRUE, or integer. +#' @param R_finite Controls variance inflation to account for finite reference LD. +#' Passed to \code{susieR::susie_rss()}. #' @param ... Additional parameters passed to susie_rss (e.g., var_y). #' @return A list with post-processed SuSiE RSS results. #' @importFrom magrittr %>% @@ -224,7 +224,7 @@ susie_rss_pipeline <- function(sumstats, LD_mat = NULL, X_mat = NULL, n = NULL, secondary_coverage = c(0.7, 0.5), signal_cutoff = 0.1, min_abs_corr = 0.8, - stochastic_ld_sample = NULL, ...) { + R_finite = NULL, ...) { analysis_method <- match.arg(analysis_method) if (!is.null(sumstats$z)) { @@ -237,7 +237,7 @@ susie_rss_pipeline <- function(sumstats, LD_mat = NULL, X_mat = NULL, n = NULL, # Common args for susie_rss_wrapper common <- list(z = z, n = n, coverage = coverage, - stochastic_ld_sample = stochastic_ld_sample, ...) + R_finite = R_finite, ...) if (!is.null(X_mat)) common$X <- X_mat else common$R <- LD_mat if (analysis_method == "single_effect") { @@ -351,7 +351,7 @@ susie_post_processor <- function(susie_output, data_x, data_y, X_scalar, y_scala if (analysis_script != "") res$analysis_script <- analysis_script if (!is.null(other_quantities)) res$other_quantities <- other_quantities if (mode == "mvsusie") { - res$context_names <- susie_output$condition_names + res$context_names <- susie_output$outcome_names } if (!is.null(data_y)) { # Mode-specific processing @@ -447,10 +447,10 @@ susie_post_processor <- function(susie_output, data_x, data_y, X_scalar, y_scala res$susie_result_trimmed$mu2 <- susie_output$mu2[eff_idx, , drop = FALSE] } if (mode == "mvsusie") { - # res$susie_result_trimmed$b1 = susie_output$b1[eff_idx, , , drop = FALSE] - # res$susie_result_trimmed$b2 = susie_output$b2[eff_idx, , , drop = FALSE] - res$susie_result_trimmed$b1_rescaled <- susie_output$b1_rescaled[eff_idx, , , drop = FALSE] - res$susie_result_trimmed$coef <- susie_output$coef + res$susie_result_trimmed$mu <- susie_output$mu[eff_idx, , , drop = FALSE] + res$susie_result_trimmed$mu2_diag <- susie_output$mu2_diag[eff_idx, , , drop = FALSE] + res$susie_result_trimmed$X_column_scale_factors <- susie_output$X_column_scale_factors + res$susie_result_trimmed$coef <- mvsusieR::coef.mvsusie(susie_output)[-1, , drop = FALSE] res$susie_result_trimmed$clfsr <- susie_output$conditional_lfsr[eff_idx, , , drop = FALSE] # other lfsr can be computed: # se_lfsr <- mvsusie_single_effect_lfsr(clfsr, alpha) diff --git a/R/twas_weights.R b/R/twas_weights.R index 415b3373..5166becb 100644 --- a/R/twas_weights.R +++ b/R/twas_weights.R @@ -634,7 +634,7 @@ twas_multivariate_weights_pipeline <- function( residual_variance = mnm_fit$mrmash_fitted$V, L = max_L, max_iter = mvsusie_max_iter, - verbosity = verbose + verbose = verbose ) ) diff --git a/R/univariate_pipeline.R b/R/univariate_pipeline.R index 63029909..13ac1133 100644 --- a/R/univariate_pipeline.R +++ b/R/univariate_pipeline.R @@ -191,7 +191,8 @@ load_study_LD <- function(ld_path, region) { #' @param impute Whether to impute missing variants via RAISS (default TRUE). #' @param impute_opts List of imputation options (rcond, R2_threshold, minimum_ld, lamb). #' @param pip_cutoff_to_skip PIP threshold for early stopping (default 0, no skip). -#' @param stochastic_ld_sample Passed to susie_rss. NULL (default), TRUE, or integer. +#' @param R_finite Controls variance inflation to account for finite reference LD. +#' Passed to \code{susieR::susie_rss()}. #' @param keep_indel Whether to keep indel variants (default TRUE). #' @param comment_string Comment character for sumstat file (default "#"). #' @param diagnostics Whether to include diagnostic info (default FALSE). @@ -211,7 +212,7 @@ rss_analysis_pipeline <- function( min_abs_corr = 0.8 ), impute = TRUE, impute_opts = list(rcond = 0.01, R2_threshold = 0.6, minimum_ld = 5, lamb = 0.01), - pip_cutoff_to_skip = 0, stochastic_ld_sample = NULL, + pip_cutoff_to_skip = 0, R_finite = NULL, keep_indel = TRUE, comment_string = "#", diagnostics = FALSE) { # Detect genotype input: single X matrix or list of X matrices (mixture panel). # susie_rss accepts X=list(X1, X2, ...) for multi-panel mixture. @@ -306,7 +307,7 @@ rss_analysis_pipeline <- function( secondary_coverage = sec_coverage, signal_cutoff = finemapping_opts$signal_cutoff, min_abs_corr = finemapping_opts$min_abs_corr, - stochastic_ld_sample = stochastic_ld_sample + R_finite = R_finite ) if (!is.null(qc_method)) { res$outlier_number <- qc_results$outlier_number diff --git a/man/rss_analysis_pipeline.Rd b/man/rss_analysis_pipeline.Rd index d1750f61..fccf216f 100644 --- a/man/rss_analysis_pipeline.Rd +++ b/man/rss_analysis_pipeline.Rd @@ -22,7 +22,7 @@ rss_analysis_pipeline( impute = TRUE, impute_opts = list(rcond = 0.01, R2_threshold = 0.6, minimum_ld = 5, lamb = 0.01), pip_cutoff_to_skip = 0, - stochastic_ld_sample = NULL, + R_finite = NULL, keep_indel = TRUE, comment_string = "#", diagnostics = FALSE @@ -65,7 +65,8 @@ signal_cutoff, min_abs_corr).} \item{pip_cutoff_to_skip}{PIP threshold for early stopping (default 0, no skip).} -\item{stochastic_ld_sample}{Passed to susie_rss. NULL (default), TRUE, or integer.} +\item{R_finite}{Controls variance inflation to account for finite reference LD. +Passed to \code{susieR::susie_rss()}.} \item{keep_indel}{Whether to keep indel variants (default TRUE).} diff --git a/man/susie_rss_pipeline.Rd b/man/susie_rss_pipeline.Rd index ead21333..2a62a8f0 100644 --- a/man/susie_rss_pipeline.Rd +++ b/man/susie_rss_pipeline.Rd @@ -17,7 +17,7 @@ susie_rss_pipeline( secondary_coverage = c(0.7, 0.5), signal_cutoff = 0.1, min_abs_corr = 0.8, - stochastic_ld_sample = NULL, + R_finite = NULL, ... ) } @@ -46,7 +46,8 @@ susie_rss_pipeline( \item{min_abs_corr}{Minimum absolute correlation for CS purity (default: 0.8).} -\item{stochastic_ld_sample}{Passed to susie_rss. NULL, TRUE, or integer.} +\item{R_finite}{Controls variance inflation to account for finite reference LD. +Passed to \code{susieR::susie_rss()}.} \item{...}{Additional parameters passed to susie_rss (e.g., var_y).} } diff --git a/man/susie_rss_wrapper.Rd b/man/susie_rss_wrapper.Rd index 65328a5f..8d1d4993 100644 --- a/man/susie_rss_wrapper.Rd +++ b/man/susie_rss_wrapper.Rd @@ -13,7 +13,7 @@ susie_rss_wrapper( max_L = 30, l_step = 5, coverage = 0.95, - stochastic_ld_sample = NULL, + R_finite = NULL, ... ) } @@ -33,9 +33,8 @@ When provided, susie_rss uses the low-rank X interface.} \item{l_step}{Step size for increasing L when the limit is reached.} -\item{stochastic_ld_sample}{Stochastic LD parameter passed to susie_rss. -NULL (default): no variance inflation. TRUE: infer sketch size from X -(requires X, not R). Integer: explicit sketch size (for R only).} +\item{R_finite}{Controls variance inflation to account for finite reference LD. +Passed to \code{susieR::susie_rss()}.} \item{...}{Extra parameters passed to susie_rss (e.g., var_y, coverage).} }