Merge pull request #183 from MangiolaLaboratory/improve-boxplot

Improve boxplot

Author: stemangiola

R/methods.R

@@ -1877,7 +1877,8 @@ sccomp_boxplot = function(
 
       pivot_wider(names_from = parameter, values_from = c(contains("c_"), contains("v_"))) |>
       unnest(count_data) |>
-      with_groups(!!.sample, ~ mutate(.x, proportion = (!!.count)/sum(!!.count)) ) 
+      with_groups(!!.sample, ~ mutate(.x, proportion = (!!.count)/sum(!!.count)) ) |> 
+      mutate(is_zero = proportion==0) 
 
   if(remove_unwanted_effects){
     .data_adjusted = 
@@ -1890,7 +1891,8 @@ sccomp_boxplot = function(
       select(-proportion) |> 
       left_join(.data_adjusted, by = join_by(!!.cell_group, !!.sample))
   }
-   
+  else 
+    message( "sccomp says: When visualising proportions, especially for complex models, consider setting `remove_unwanted_effects=TRUE`. This will adjust the proportions, preserving only the observed effect.")
 
   # If I don't have outliers add them
   if(!"outlier" %in% colnames(data_proportion)) data_proportion = data_proportion |> mutate(outlier = FALSE) 
@@ -2013,13 +2015,9 @@ else {
         # If discrete
         else 
           my_plot = 
-            plot_boxplot(
+            sccomp_boxplot(
               .data = x,
-              data_proportion = data_proportion,
-              factor_of_interest = .x,
-              .cell_group = !!.cell_group,
-              .sample =  !!.sample,
-              my_theme = multipanel_theme,
+              factor = .x,
               significance_threshold = significance_threshold
             ) 
 

R/utilities.R

@@ -2028,7 +2028,10 @@ plot_2D_intervals = function(
 #' @noRd
 plot_boxplot = function(
     .data, data_proportion, factor_of_interest, .cell_group,
-    .sample, significance_threshold = 0.05, my_theme, remove_unwanted_effects = FALSE
+    .sample, 
+    significance_threshold = 0.05, 
+    my_theme, 
+    remove_unwanted_effects = FALSE
 ){
 
   # Define the variables as NULL to avoid CRAN NOTES
@@ -2123,15 +2126,15 @@ plot_boxplot = function(
       simulated_proportion =
         .data |>
         sccomp_replicate(formula_composition = formula_composition, number_of_draws = 100) |>
-        left_join(data_proportion %>% distinct(!!as.symbol(factor_of_interest), !!.sample, !!.cell_group))
+        left_join(data_proportion %>% distinct(!!as.symbol(factor_of_interest), !!.sample, !!.cell_group, is_zero))
 
     my_boxplot = my_boxplot +
 
       # Add boxplot for simulated proportions
       stat_summary(
         aes(!!as.symbol(factor_of_interest), (generated_proportions)),
         fun.data = calc_boxplot_stat, geom="boxplot",
-        outlier.shape = NA, outlier.color = NA,outlier.size = 0,
+        outlier.shape = NA, outlier.color = NA, outlier.size = 0,
         fatten = 0.5, lwd=0.2,
         data =
           simulated_proportion %>%
@@ -2180,7 +2183,7 @@ plot_boxplot = function(
 
     # Add jittered points for individual data
     geom_jitter(
-      aes(!!as.symbol(factor_of_interest), proportion, shape=outlier, color=outlier,  group=!!as.symbol(factor_of_interest)),
+      aes(!!as.symbol(factor_of_interest), proportion, shape=is_zero, color=outlier,  group=!!as.symbol(factor_of_interest)),
       data = data_proportion,
       position=position_jitterdodge(jitter.height = 0, jitter.width = 0.2),
       size = 0.5
@@ -2193,13 +2196,13 @@ plot_boxplot = function(
       nrow = 4
     ) +
     scale_color_manual(values = c("black", "#e11f28")) +
+    scale_shape_manual(values = c(16, 21)) +
     scale_y_continuous(trans=S_sqrt_trans(), labels = dropLeadingZero) +
     scale_fill_discrete(na.value = "white") +
     xlab("Biological condition") +
     ylab("Cell-group proportion") +
-    guides(color="none", alpha="none", size="none") +
+    guides( alpha="none", size="none") +
     labs(fill="Significant difference") +
-    ggtitle("Note: Be careful judging significance (or outliers) visually for lowly abundant cell groups. \nVisualising proportion hides the uncertainty characteristic of count data, that a count-based statistical model can estimate.") +
     my_theme +
     theme(axis.text.x =  element_text(angle=20, hjust = 1), title = element_text(size = 3))
 }

README.md

@@ -170,6 +170,13 @@ sccomp_result |>
   sccomp_boxplot(factor = "type")
 ```
 
+You can plot proportions adjusted for unwanted effects. This is helpful especially for complex models, where multiple factors can significantly impact the proportions.
+
+```{r, eval = instantiate::stan_cmdstan_exists()}
+sccomp_result |> 
+  sccomp_boxplot(factor = "type", remove_unwanted_effects = TRUE)
+```
+
 A plot of estimates of differential composition (c_) on the x-axis and differential variability (v_) on the y-axis. The error bars represent 95% credible intervals. The dashed lines represent the minimal effect that the hypothesis test is based on. An effect is labelled as significant if it exceeds the minimal effect according to the 95% credible interval. Facets represent the covariates in the model.
 
 ```{r, eval = instantiate::stan_cmdstan_exists()}

inst/figures/unnamed-chunk-12-1.png

@@ -219,6 +219,13 @@ sccomp_result |>
   sccomp_boxplot(factor = "type")
 ```
 
+You can plot proportions adjusted for unwanted effects. This is helpful especially for complex models, where multiple factors can significantly impact the proportions.
+
+```{r, eval = instantiate::stan_cmdstan_exists()}
+sccomp_result |> 
+  sccomp_boxplot(factor = "type", remove_unwanted_effects = TRUE)
+```
+
 A plot of estimates of differential composition (c_) on the x-axis and differential variability (v_) on the y-axis. The error bars represent 95% credible intervals. The dashed lines represent the minimal effect that the hypothesis test is based on. An effect is labelled as significant if it exceeds the minimal effect according to the 95% credible interval. Facets represent the covariates in the model.
 
 ```{r, eval = instantiate::stan_cmdstan_exists()}