Better plots and vignette

- annotation of paths in the tree
- growth model vignette clarified and issue opened

R/annotate_forest.R

@@ -0,0 +1,93 @@
+#' Annotate a plot of cell divisions.
+#'
+#' @description
+#' It annotates a plot of cell divisions with information from sampling
+#' times and MRCAs for all available samples
+#'
+#' @param tree_plot The output of `plot_forest`.
+#' @param forest The original forest object from which the input to
+#'        `plot_forest`
+#' has been derived.
+#' @param samples If `TRUE` it annotates samples.
+#' @param MRCAs If `TRUE` it annotates MRCAs
+#'
+#' @return A `ggraph` tree plot.
+#' @export
+#'
+#' @examples
+#' sim <- new(Simulation)
+#' sim$add_genotype(name = "A", growth_rates = 0.08, death_rates = 0.01)
+#' sim$place_cell("A", 500, 500)
+#' sim$run_up_to_time(60)
+#' sim$sample_cells("MySample", c(500, 500), c(510, 510))
+#' forest = sim$get_samples_forest()
+#' forest$get_samples_info()
+#' tree_plot = plot_forest(forest)
+#' annotate_forest(tree_plot, forest)
+annotate_forest <- function(tree_plot, forest, samples = TRUE, MRCAs = TRUE) {
+  samples_info <- forest$get_samples_info()
+
+  # Sampling times
+  if (samples) {
+
+    max_Y <- max(tree_plot$data$y, na.rm = TRUE)
+
+    tree_plot <- tree_plot +
+      ggplot2::geom_hline(
+        yintercept = max_Y - samples_info$time,
+        color = "indianred3",
+        linetype = "dashed",
+        linewidth = .3
+      )
+  }
+
+  # MRCAs
+  if(MRCAs) {
+    sample_names <- samples_info %>% pull(.data$name)
+
+    MRCAs_cells <- lapply(sample_names,
+                          function(s){
+                            forest$get_coalescent_cells(
+                              forest$get_nodes() %>%
+                                dplyr::filter(sample %in% s) %>%
+                                dplyr::pull(.data$cell_id)
+                            ) %>%
+                              dplyr::mutate(sample = s)
+                          }) %>%
+      Reduce(f = dplyr::bind_rows) %>%
+      dplyr::group_by(.data$cell_id) %>%
+      dplyr::mutate(
+        cell_id = paste(.data$cell_id)
+      ) %>%
+      dplyr::summarise(
+        label = paste0("    ", .data$sample, collapse = "\n")
+      )
+
+    layout <- tree_plot$data %>%
+      dplyr::select(.data$x, .data$y, .data$name) %>%
+      dplyr::mutate(cell_id = paste(.data$name)) %>%
+      dplyr::filter(.data$name %in% MRCAs_cells$cell_id) %>%
+      dplyr::left_join(MRCAs_cells, by = "cell_id")
+
+    tree_plot <-
+      tree_plot +
+      ggplot2::geom_point(
+        data = layout,
+        ggplot2::aes(x = .data$x, y = .data$y),
+        color = "purple3",
+        size = 3,
+        pch = 21
+      ) +
+      ggplot2::geom_text(
+        data = layout,
+        ggplot2::aes(x = .data$x, y = .data$y,
+                     label = .data$label),
+        color = "purple3",
+        size = 3,
+        hjust = 0,
+        vjust = 1
+      )
+  }
+
+  tree_plot
+}

R/plot_forest.R

@@ -22,6 +22,8 @@
 #' set of trees connected to a generic wildtype cell.
 #'
 #' @param forest The samples forest to be plot.
+#' @param highlight_sample If a sample name, the path from root to the sampled
+#' cells in the sample is highlighted. If `NULL` (default), nothing is highlighted.
 #'
 #' @return A `ggraph` tree plot.
 #' @export
@@ -34,7 +36,7 @@
 #' sim$sample_cells("MySample", c(500, 500), c(510, 510))
 #' forest <- sim$get_samples_forest()
 #' plot_forest(forest)
-plot_forest <- function(forest) {
+plot_forest <- function(forest, highlight_sample = NULL) {
   stopifnot(inherits(forest, "Rcpp_SamplesForest"))
 
   nodes <- forest$get_nodes()
@@ -54,11 +56,19 @@ plot_forest <- function(forest) {
       dplyr::mutate(
         from = ifelse(is.na(.data$from), "WT", .data$from),
         species = paste0(.data$genotype, .data$epistate),
-        sample = ifelse(is.na(.data$sample), "N/A", .data$sample)
+        sample = ifelse(is.na(.data$sample), "N/A", .data$sample),
+        highlight = FALSE
       )
+
+    # Highlight is optional
+    if(!is.null(highlight_sample))
+    {
+      highlight = paths_to_sample(forest_data, highlight_sample)
+      forest_data$highlight = forest_data$to %in% highlight$to
+    }
 
     edges <- forest_data %>%
-      dplyr::select(.data$from, .data$to)
+      dplyr::select("from", "to", "highlight")
 
     # Create a tidygraph object
     graph <- tidygraph::as_tbl_graph(edges, directed = TRUE)
@@ -100,15 +110,31 @@ plot_forest <- function(forest) {
     point_size <- c(.5, rep(1, nsamples))
     names(point_size) <- c("N/A", forest$get_samples_info() %>%
                             pull(.data$name))
-
-    # Plot the graph with color-coded nodes
-    ggraph::ggraph(layout, "tree") +
-      ggraph::geom_edge_link(edge_width = .1) +
+
+    # highlight_edges = paths_to_sample(layout, "S_2_2")
+
+    # layout$highlight = layout$name %in% highlight_edges$name
+
+    # Plot the graph
+    graph_plot = ggraph::ggraph(layout, "tree") +
+      ggraph::geom_edge_link(edge_width = .1, 
+                             aes(edge_color = ifelse(highlight, 'indianred3', 'black'))
+                             )
+
+    graph_plot +
       ggraph::geom_node_point(ggplot2::aes(color = .data$species,
                                            shape = ifelse(is.na(.data$sample),
                                                           "N/A",
                                                           .data$sample),
-                                           size = .data$sample)) +
+                                           size = .data$sample)
+                              ) +
+      # ggraph::geom_edge_link(edge_width = .1, ggplot2::aes(color = birth_time < 46)) +
+      # ggraph::geom_edge_link(data = highlight_edges, color = 'red') +
+      # ggraph::geom_node_point(ggplot2::aes(color = .data$species,
+      #                                      shape = ifelse(is.na(.data$sample),
+      #                                                     "N/A",
+      #                                                     .data$sample),
+      #                                      size = .data$sample)) 
       ggplot2::scale_color_manual(values = species_colors) +
       # scale_shape_manual(values = c("+" = 16, "-" = 17)) +
       # theme_void() +
@@ -140,96 +166,27 @@ plot_forest <- function(forest) {
   }
 }
 
-#' Annotate a plot of cell divisions.
-#'
-#' @description
-#' It annotates a plot of cell divisions with information from sampling
-#' times and MRCAs for all available samples
-#'
-#' @param tree_plot The output of `plot_forest`.
-#' @param forest The original forest object from which the input to
-#'        `plot_forest`
-#' has been derived.
-#' @param samples If `TRUE` it annotates samples.
-#' @param MRCAs If `TRUE` it annotates MRCAs
-#'
-#' @return A `ggraph` tree plot.
-#' @export
-#'
-#' @examples
-#' sim <- new(Simulation)
-#' sim$add_genotype(name = "A", growth_rates = 0.08, death_rates = 0.01)
-#' sim$place_cell("A", 500, 500)
-#' sim$run_up_to_time(60)
-#' sim$sample_cells("MySample", c(500, 500), c(510, 510))
-#' forest = sim$get_samples_forest()
-#' forest$get_samples_info()
-#' tree_plot = plot_forest(forest)
-#' annotate_forest(tree_plot, forest)
-annotate_forest <- function(tree_plot, forest, samples = TRUE, MRCAs = TRUE) {
-  samples_info <- forest$get_samples_info()
-
-  # Sampling times
-  if (samples) {
-
-    max_Y <- max(tree_plot$data$y, na.rm = TRUE)
-
-    tree_plot <- tree_plot +
-      ggplot2::geom_hline(
-        yintercept = max_Y - samples_info$time,
-        color = "indianred3",
-        linetype = "dashed",
-        linewidth = .3
-      )
-  }
-
-  # MRCAs
-  if(MRCAs) {
-    sample_names <- samples_info %>% pull(.data$name)
-
-    MRCAs_cells <- lapply(sample_names,
-                          function(s){
-                            forest$get_coalescent_cells(
-                              forest$get_nodes() %>%
-                                dplyr::filter(sample %in% s) %>%
-                                dplyr::pull(.data$cell_id)
-                            ) %>%
-                              dplyr::mutate(sample = s)
-                          }) %>%
-      Reduce(f = dplyr::bind_rows) %>%
-      dplyr::group_by(.data$cell_id) %>%
-      dplyr::mutate(
-        cell_id = paste(.data$cell_id)
-      ) %>%
-      dplyr::summarise(
-        label = paste0("    ", .data$sample, collapse = "\n")
-      )
 
-    layout <- tree_plot$data %>%
-      dplyr::select(.data$x, .data$y, .data$name) %>%
-      dplyr::mutate(cell_id = paste(.data$name)) %>%
-      dplyr::filter(.data$name %in% MRCAs_cells$cell_id) %>%
-      dplyr::left_join(MRCAs_cells, by = "cell_id")
-
-    tree_plot <-
-      tree_plot +
-      ggplot2::geom_point(
-        data = layout,
-        ggplot2::aes(x = .data$x, y = .data$y),
-        color = "purple3",
-        size = 3,
-        pch = 21
-      ) +
-      ggplot2::geom_text(
-        data = layout,
-        ggplot2::aes(x = .data$x, y = .data$y,
-                     label = .data$label),
-        color = "purple3",
-        size = 3,
-        hjust = 0,
-        vjust = 1
-      )
+paths_to_sample = function(forest_data, sample){
+
+  to = forest_data %>% dplyr::filter(sample == !!sample) %>% pull(to)
+
+  Q = NULL
+
+  while(length(to) > 0)
+  {
+    # New element to inspect
+    to_head = to[1]
+    to = to[-1]
+
+    # Forward star 
+    to_add = forest_data %>%  dplyr::filter(to == to_head)
+    Q = dplyr::bind_rows(Q, to_add) %>% dplyr::distinct()
+
+    # Recursion
+    to = c(to, to_add %>% dplyr::pull(from)) %>% unique()
   }
-
-  tree_plot
+  
+  Q
 }
+

vignettes/a03_growth_models.Rmd

@@ -21,7 +21,7 @@ knitr::opts_chunk$set(
 library(rRACES)
 ```
 
-rRACES/RACS supports two growth models: the "border"-growth model and the homogeneous-growth
+rRACES/RACES supports two growth models: the "border"-growth model and the homogeneous-growth
 model. 
 
 The former exclusively admits duplication of cells that have space to duplicate 
@@ -116,7 +116,7 @@ plot_muller(sim)
 
 Finally, let us build the ancestor forest of the samples.
 
-```{r}
+```{r, fig.width=14}
 library(dplyr)
 library(ggplot2)
 
@@ -126,6 +126,24 @@ plot_forest(forest) %>%
   annotate_forest(forest)
 ```
 
+We use a special `highlight` parameter to add the edges connecting cells in each sample.
+```{r, fig.width=14}
+plot_forest(forest, highlight ='S_1_1') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_1_2') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_2_1') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_2_2') %>%
+  annotate_forest(forest)
+```
+
+
+
+
 ### "Border" Growth
 
 We now build a simulation and switch to the "border" growth model.
@@ -199,16 +217,28 @@ plot_muller(sim)
 
 Once more, let us build the ancestor forest of the samples.
 
-```{r}
-library(dplyr)
-library(ggplot2)
-
+```{r, fig.width=14}
 forest <- sim$get_samples_forest()
 
 plot_forest(forest) %>%
   annotate_forest(forest)
 ```
 
+We use a special `highlight` parameter to add the edges connecting cells in each sample.
+```{r, fig.width=14}
+plot_forest(forest, highlight ='S_1_1') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_1_2') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_2_1') %>%
+  annotate_forest(forest)
+
+plot_forest(forest, highlight ='S_2_2') %>%
+  annotate_forest(forest)
+```
+
 It is easy to see the differences in the analoguous plots of the above examples.
 
 First of all, the cancer growth is slower when subject to the "border"-growth model