diff --git a/R/functions_multi_beta_binomial.R b/R/functions_multi_beta_binomial.R index 1fef5cee..4c379bc8 100644 --- a/R/functions_multi_beta_binomial.R +++ b/R/functions_multi_beta_binomial.R @@ -496,27 +496,14 @@ multi_beta_binomial_glm = function(.data, add_attr(.cell_group, ".cell_group") |> add_attr(.count, ".count") |> add_attr(check_outliers, "check_outliers") |> - + add_attr(formula_composition, "formula_composition") |> + add_attr(formula_variability, "formula_variability") |> test_contrasts( contrasts = contrasts, percent_false_positive = percent_false_positive, test_composition_above_logit_fold_change = test_composition_above_logit_fold_change - ) |> - - # Add osme attributes - add_attr(get_mean_precision_association(estimates_list$fit), "mean_concentration_association") %>% - when(pass_fit ~ add_attr(., estimates_list$fit, "fit"), ~ (.)) |> - - # Attach association mean concentration - add_attr(estimates_list$data_for_model, "model_input") |> - add_attr(estimates_list$truncation_df2, "truncation_df2") |> - add_attr(.sample, ".sample") |> - add_attr(.cell_group, ".cell_group") |> - add_attr(.count, ".count") |> - add_attr(check_outliers, "check_outliers") |> - add_attr(formula_composition, "formula_composition") |> - add_attr(formula_variability, "formula_variability") + ) } diff --git a/R/methods.R b/R/methods.R index 34608a1e..8f6d4738 100644 --- a/R/methods.R +++ b/R/methods.R @@ -594,6 +594,7 @@ sccomp_glm_data_frame_counts = function(.data, #' @param contrasts A vector of character strings. For example if your formula is `~ 0 + treatment` and the factor treatment has values `yes` and `no`, your contrast could be "constrasts = c(treatmentyes - treatmentno)". #' @param percent_false_positive A real between 0 and 100 non included. This used to identify outliers with a specific false positive rate. #' @param test_composition_above_logit_fold_change A positive integer. It is the effect threshold used for the hypothesis test. A value of 0.2 correspond to a change in cell proportion of 10% for a cell type with baseline proportion of 50%. That is, a cell type goes from 45% to 50%. When the baseline proportion is closer to 0 or 1 this effect thrshold has consistent value in the logit uncontrained scale. +#' @param pass_fit A boolean. Whether to pass the Stan fit as attribute in the output. Because the Stan fit can be very large, setting this to FALSE can be used to lower the memory imprint to save the output. #' #' @return A nested tibble `tbl` with cell_group-wise statistics #' @@ -616,7 +617,8 @@ sccomp_glm_data_frame_counts = function(.data, test_contrasts <- function(.data, contrasts = NULL, percent_false_positive = 5, - test_composition_above_logit_fold_change = 0.2) { + test_composition_above_logit_fold_change = 0.2, + pass_fit = TRUE) { UseMethod("test_contrasts", .data) } @@ -626,7 +628,8 @@ test_contrasts <- function(.data, test_contrasts.data.frame = function(.data, contrasts = NULL, percent_false_positive = 5, - test_composition_above_logit_fold_change = 0.2){ + test_composition_above_logit_fold_change = 0.2, + pass_fit = TRUE){ .sample = .data |> attr(".sample") @@ -640,7 +643,7 @@ test_contrasts.data.frame = function(.data, abundance_CI = get_abundance_contrast_draws(.data, contrasts) %>% - # If my constrasts do not match my model. I have to do something more elegant. + # If my contrasts do not match my model. I have to do something more elegant. when( "parameter" %in% colnames(.) ~ draws_to_statistics( ., @@ -702,7 +705,21 @@ test_contrasts.data.frame = function(.data, nest(count_data = -!!.cell_group), by = quo_name(.cell_group) ) - ) + ) |> + + # Add back attributes + add_attr(.data |> attr("fit") |> get_mean_precision_association(), "mean_concentration_association") %>% + when(pass_fit ~ add_attr(., .data |> attr("fit") , "fit"), ~ (.)) |> + + # Attach association mean concentration + add_attr(.data |> attr("model_input") , "model_input") |> + add_attr(.data |> attr("truncation_df2"), "truncation_df2") |> + add_attr(.sample, ".sample") |> + add_attr(.cell_group, ".cell_group") |> + add_attr(.count, ".count") |> + add_attr(check_outliers, "check_outliers") |> + add_attr(.data |> attr("formula_composition"), "formula_composition") |> + add_attr(.data |> attr("formula_variability"), "formula_variability") } #' sccomp_replicate diff --git a/man/test_contrasts.Rd b/man/test_contrasts.Rd index b61c9b83..2c9f00a3 100644 --- a/man/test_contrasts.Rd +++ b/man/test_contrasts.Rd @@ -8,7 +8,8 @@ test_contrasts( .data, contrasts = NULL, percent_false_positive = 5, - test_composition_above_logit_fold_change = 0.2 + test_composition_above_logit_fold_change = 0.2, + pass_fit = TRUE ) } \arguments{ @@ -19,6 +20,8 @@ test_contrasts( \item{percent_false_positive}{A real between 0 and 100 non included. This used to identify outliers with a specific false positive rate.} \item{test_composition_above_logit_fold_change}{A positive integer. It is the effect threshold used for the hypothesis test. A value of 0.2 correspond to a change in cell proportion of 10\% for a cell type with baseline proportion of 50\%. That is, a cell type goes from 45\% to 50\%. When the baseline proportion is closer to 0 or 1 this effect thrshold has consistent value in the logit uncontrained scale.} + +\item{pass_fit}{A boolean. Whether to pass the Stan fit as attribute in the output. Because the Stan fit can be very large, setting this to FALSE can be used to lower the memory imprint to save the output.} } \value{ A nested tibble \code{tbl} with cell_group-wise statistics diff --git a/tests/testthat/test-sccomp_.R b/tests/testthat/test-sccomp_.R index c09a19ae..19569c74 100644 --- a/tests/testthat/test-sccomp_.R +++ b/tests/testthat/test-sccomp_.R @@ -427,3 +427,25 @@ test_that("plot test variability",{ }) + +test_that("test constrasts",{ + + + estimate = + seurat_obj |> + sccomp_glm( + formula_composition = ~ type , + formula_variability = ~ 1, + sample, cell_group, + check_outliers = FALSE, + approximate_posterior_inference = FALSE, + cores = 1, + mcmc_seed = 42 + ) + + new_test = + estimate |> + test_contrasts() |> + test_contrasts() + +})