avg coverage

Author: showteeth

NEWS.md

@@ -6,7 +6,8 @@
 
 ## Minor changes
 * `geom_coverage` and `ggcoverage` supports changing shared type of y-axis scales across facets (with `facet.y.scale` parameter).
-* `geom_coverage` and `ggcoverage` supports joint view of multiple tracks.
+* `geom_coverage` and `ggcoverage` supports joint view of multiple tracks (with `plot.type` parameter).
+* `geom_coverage` and `ggcoverage` supports joint view of average coverage of multiple tracks (with `joint.avg` parameter).
 * `LoadTrackFile` loads bam directly instead of converting bam to bigwig and loading the bigwig.
 
 -------------

R/geom_coverage.R

@@ -7,7 +7,8 @@
 #' @param single.nuc Logical value, whether to visualize at single nucleotide level (use bar plot). Default: FALSE.
 #' @param plot.type The type of the plot, choose from facet (separate plot for every sample) and
 #' joint (combine all sample in a single plot). Default: facet.
-#' @param facet.key Sample type key to create coverage plot. Default: Type.
+#' @param facet.key Sample type key to create coverage plot. Used in both facet and joint plot. Default: Type.
+#' @param joint.avg Logical value, whether to show average coverage across \code{group.key}. Default: FALSE.
 #' @param facet.order The order of coverage plot. Default: NULL.
 #' @param facet.color The color of sample text. Default: NULL (select automatically).
 #' @param facet.y.scale The shared type of y-axis scales across facets, choose from free (facets have different y-axis scales),
@@ -52,7 +53,7 @@
 #' #   geom_coverage(data = track.df, color = "auto", mark.region = NULL)
 geom_coverage <- function(data, mapping = NULL, color = NULL, rect.color = NA,
                           single.nuc = FALSE, plot.type = c("facet", "joint"),
-                          facet.key = "Type", facet.order = NULL,
+                          facet.key = "Type", joint.avg = FALSE, facet.order = NULL,
                           facet.color = NULL, facet.y.scale = c("free", "fixed"),
                           group.key = "Group", range.size = 3, range.position = c("in", "out"),
                           mark.region = NULL, mark.color = "grey", mark.alpha = 0.5,
@@ -231,28 +232,60 @@ geom_coverage <- function(data, mapping = NULL, color = NULL, rect.color = NA,
       plot.ele <- append(plot.ele, region.range)
     }
   } else if (plot.type == "joint") {
-    # preocess data
     data.split <- split(x = data, f = data[, facet.key])
-    data.split.final <- lapply(data.split, function(ds) {
-      ds <- ds[order(ds$start), ]
-      ds.ll <- ds[nrow(ds), ]
-      if (ds.ll$start < ds.ll$end) {
-        ds.ll.new <- ds.ll
-        ds.ll.new$start <- ds.ll.new$end
-        ds.final <- rbind(ds, ds.ll.new) %>% as.data.frame()
-      } else {
-        ds.final <- ds
-      }
-      rownames(ds.final) <- NULL
-      return(ds.final)
-    })
-    # get final data
-    data.final <- do.call(rbind, data.split.final) %>% as.data.frame()
-    rownames(data.final) <- NULL
+    if (joint.avg) {
+      # process data
+      data.list <- lapply(data.split, function(ds) {
+        ds.list <- apply(ds, 1, function(x) {
+          x.start <- seq(x["start"], x["end"])
+          x.len <- length(x.start)
+          x.df <- data.frame(
+            seqnames = rep(x["seqnames"], x.len), start = x.start,
+            score = rep(x["score"], x.len)
+          )
+          x.df[, facet.key] <- x[facet.key]
+          x.df[, group.key] <- x[group.key]
+          x.df
+        })
+        ds.final <- do.call(rbind, ds.list)
+        rownames(ds.final) <- NULL
+        return(ds.final)
+      })
+      data.final <- do.call(rbind, data.list)
+      rownames(data.final) <- NULL
+      # as numeric
+      data.final$score <- as.numeric(data.final$score)
+      data.final$start <- as.numeric(data.final$start)
+      # get average
+      data.final.mean <- data.final %>%
+        dplyr::group_by(.data[[group.key]], .data[["start"]]) %>%
+        dplyr::summarise(score = mean(.data[["score"]]))
+      # create avg plot
+      avg.plot <- geom_line(data = data.final.mean, mapping = mapping)
+      plot.ele <- list(avg.plot)
+    } else {
+      # preocess data
+      data.split.final <- lapply(data.split, function(ds) {
+        ds <- ds[order(ds$start), ]
+        ds.ll <- ds[nrow(ds), ]
+        if (ds.ll$start < ds.ll$end) {
+          ds.ll.new <- ds.ll
+          ds.ll.new$start <- ds.ll.new$end
+          ds.final <- rbind(ds, ds.ll.new) %>% as.data.frame()
+        } else {
+          ds.final <- ds
+        }
+        rownames(ds.final) <- NULL
+        return(ds.final)
+      })
+      # get final data
+      data.final <- do.call(rbind, data.split.final) %>% as.data.frame()
+      rownames(data.final) <- NULL
 
-    # create step plot
-    region.step <- geom_step(data = data.final, mapping = mapping)
-    plot.ele <- list(region.step)
+      # create step plot
+      region.step <- geom_step(data = data.final, mapping = mapping)
+      plot.ele <- list(region.step)
+    }
   }
 
   # add rect

R/ggcoverage.R

@@ -8,6 +8,7 @@
 #' @param plot.type The type of the plot, choose from facet (separate plot for every sample) and
 #' joint (combine all sample in a single plot). Default: facet.
 #' @param facet.key Sample type key to create coverage plot. Default: Type.
+#' @param joint.avg Logical value, whether to show average coverage across \code{group.key}. Default: FALSE.
 #' @param facet.order The order of Coverage plot. Default: NULL.
 #' @param facet.color The color of sample text. Default: NULL (select automatically).
 #' @param facet.y.scale The shared type of y-axis scales across facets, choose from free (facets have different y-axis scales),
@@ -30,7 +31,7 @@
 #' @importFrom scales comma
 #' @importFrom grDevices colorRampPalette
 #' @importFrom RColorBrewer brewer.pal
-#' @importFrom rlang as_label
+#' @importFrom rlang as_label .data
 #' @importFrom stats as.formula
 #' @importFrom ggh4x facet_wrap2 strip_themed
 #' @importFrom dplyr group_by summarise
@@ -53,7 +54,7 @@
 #' # gtf.gr <- rtracklayer::import.gff(con = gtf.file, format = "gtf")
 #' # ggcoverage(data = track.df, color = "auto", range.position = "out")
 ggcoverage <- function(data, single.nuc = FALSE, mapping = NULL, color = NULL,
-                       rect.color = NA, plot.type = c("facet", "joint"), facet.key = "Type",
+                       rect.color = NA, plot.type = c("facet", "joint"), facet.key = "Type", joint.avg = FALSE,
                        facet.order = NULL, facet.color = NULL, facet.y.scale = c("free", "fixed"),
                        group.key = "Group", range.size = 3, range.position = c("in", "out"), plot.space = 0.2,
                        mark.region = NULL, mark.color = "grey", mark.alpha = 0.5, show.mark.label = TRUE, mark.label.size = 4) {
@@ -65,7 +66,7 @@ ggcoverage <- function(data, single.nuc = FALSE, mapping = NULL, color = NULL,
     geom_coverage(
       data = data, mapping = mapping, color = color, rect.color = rect.color,
       single.nuc = single.nuc, plot.type = plot.type,
-      facet.key = facet.key, facet.order = facet.order,
+      facet.key = facet.key, joint.avg = joint.avg, facet.order = facet.order,
       facet.color = facet.color, facet.y.scale = facet.y.scale,
       group.key = group.key, range.size = range.size, range.position = range.position,
       mark.region = mark.region, mark.color = mark.color, mark.alpha = mark.alpha, show.mark.label = show.mark.label,

README.Rmd

@@ -119,9 +119,19 @@ The basic coverage plot has **two types**:
 * **joint**: Visualize all tracks in a single plot.
 
 #### joint view
+Create line plot for **every sample** (`facet.key = "Type"`) and color by **every sample** (`group.key = "Type"`):
 ```{r basic_coverage_joint, warning=FALSE, fig.height = 4, fig.width = 12, fig.align = "center"}
 basic.coverage = ggcoverage(data = track.df, color = "auto", 
-                            plot.type = "joint", group.key = "Type",
+                            plot.type = "joint", facet.key = "Type", group.key = "Type",
+                            mark.region = mark.region, range.position = "out")
+basic.coverage
+```
+
+Create **group average line plot** (sample is indicated by `facet.key = "Type"`, group is indicated by `group.key = "Group"`):
+```{r basic_coverage_joint_avg, warning=FALSE, fig.height = 4, fig.width = 12, fig.align = "center"}
+basic.coverage = ggcoverage(data = track.df, color = "auto", 
+                            plot.type = "joint", facet.key = "Type", group.key = "Group", 
+                            joint.avg = TRUE,
                             mark.region = mark.region, range.position = "out")
 basic.coverage
 ```

README.md

@@ -158,15 +158,31 @@ The basic coverage plot has **two types**:
 
 #### joint view
 
+Create line plot for **every sample** (`facet.key = "Type"`) and color
+by **every sample** (`group.key = "Type"`):
+
 ``` r
 basic.coverage = ggcoverage(data = track.df, color = "auto", 
-                            plot.type = "joint", group.key = "Type",
+                            plot.type = "joint", facet.key = "Type", group.key = "Type",
                             mark.region = mark.region, range.position = "out")
 basic.coverage
 ```
 
 <img src="man/figures/README-basic_coverage_joint-1.png" width="100%" style="display: block; margin: auto;" />
 
+Create **group average line plot** (sample is indicated by
+`facet.key = "Type"`, group is indicated by `group.key = "Group"`):
+
+``` r
+basic.coverage = ggcoverage(data = track.df, color = "auto", 
+                            plot.type = "joint", facet.key = "Type", group.key = "Group", 
+                            joint.avg = TRUE,
+                            mark.region = mark.region, range.position = "out")
+basic.coverage
+```
+
+<img src="man/figures/README-basic_coverage_joint_avg-1.png" width="100%" style="display: block; margin: auto;" />
+
 #### facet view
 
 ``` r

man/figures/README-basic_coverage_joint_avg-1.png

@@ -136,9 +136,10 @@ The basic coverage plot has **two types**:
 * **joint**: Visualize all tracks in a single plot.
 
 #### joint view
+Create line plot for **every sample** (`facet.key = "Type"`) and color by **every sample** (`group.key = "Type"`):
 ```{r basic_coverage_joint, eval=FALSE}
 basic.coverage = ggcoverage(data = track.df, color = "auto", 
-                            plot.type = "joint", group.key = "Type",
+                            plot.type = "joint", facet.key = "Type", group.key = "Type",
                             mark.region = mark.region, range.position = "out")
 basic.coverage
 ```
@@ -147,6 +148,19 @@ basic.coverage
 knitr::include_graphics("../man/figures/README-basic_coverage_joint-1.png")
 ```
 
+Create **group average line plot** (sample is indicated by `facet.key = "Type"`, group is indicated by `group.key = "Group"`):
+```{r basic_coverage_joint_avg, eval=FALSE}
+basic.coverage = ggcoverage(data = track.df, color = "auto", 
+                            plot.type = "joint", facet.key = "Type", group.key = "Group", 
+                            joint.avg = TRUE,
+                            mark.region = mark.region, range.position = "out")
+basic.coverage
+```
+
+```{r basic_coverage_joint_avg_plot, echo=FALSE, fig.height = 4, fig.width = 12, fig.align = "center"}
+knitr::include_graphics("../man/figures/README-basic_coverage_joint_avg-1.png")
+```
+
 #### facet view
 ```{r basic_coverage, eval=FALSE}
 basic.coverage = ggcoverage(data = track.df, color = "auto", plot.type = "facet",